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INTRODUCTION 

This NASA Grant (NGL 22-009-187) covers the support 
of several MIT faculty members and some of their graduate 
students in research efforts concerned with the Moon. 
Primary emphasis is on the interpretation of lunar data 
developed during the Apollo Program with additional data 
being used as needed for greater understanding. The 
Attached Annual Report (covering the period 1 March 1972 
through 28 February 1973) is presented in mixed form, 
including both finished papers and brief progress reports 
describing the efforts of various individuals being 
supported by the grant. 



X-RAY EMISSION SPECTRA AND MOLECULAR ORBITAL 
BONDING MODELS OF LUNAR MATERIALS 


by 


J.A. Tossell, D.J. Vaughan and R.G. Burns 
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Objectives 


The aim of the project is to calculate the electronic struc- 
tures of important coordination polyhedra occurring in lunar min- 
erals. These include the SiO^ tetrahedron in olivines, the 
Si 2 0^^ group in pyroxenes, the MgOg^ and FeOg^ octahedra in 

. . - 4 ' 

ferromagnesian silicates and the FeO^ tetrahedron m spinels. 

-9 -10 

We plan to extend the computations to the TiOg and CrOg octa- 
hedra (containing Ti and Cr ) believed to occur in the ferro- 
magnesian silicate minerals which crystallized under very low 
oxygen fugacities on the moon. From accurate molecular orbital 
(MO) calculations, the energies (wavelengths) at which allowed 


and forbidden crystal-field bands, as well as metal -metal and 
metal-oxygen charge transfer transitions, occur may be predicted. 
Consequently, absorption spectral assignments may be clarified 
and information obtained on the oxidation states and crystal 
chemistry of Fe,Ti and Cr in lunar and terrestrial minerals. 

A second approach has been to measure the soft x-ray emission 
spectra of Fe , Ti and Cr in lunar and terrestrial minerals. The 
positions of the peaks in the x-ray spectra not only can be 
predicted from the MO calculations and serve as a check on their 
validity, but they also provide an independent method for deter- 
mining the oxidation states and coordination numbers of the ele- 
ments in the lunar minerals. We have obtained excellent results 
for iron in terrestrial and lunar minerals, and are currently 
extending the measurements to Ti and Cr. 
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Molecular Orbital Project 
J.A. Tossell 

'■ . ' . and . 

D.J. Vaughan . - 

The molecular orbital (MO) interpretation of the x-ray emis- 
sion spectra of minerals has been placed on a firm quantitative 
basis by a series of calculations on the spectra of Mg, A1 and 
Si oxides. This work has resulted in two publications and one 
abstract, copies of which are enclosed. Our current work employs 
a new MO method for the study of oxidic minerals of both transi- 
tion and non-transition metals. Our first paper using this method 
on the spectra and electronic structure of quartz, is enclosed. 
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MOLECULAR ORBITAL INTERPRETATION OF 
X-RAY EMISSION AND ESCA SPECTRAL SHIFTS 

IN SILICATES* 

.1. a. tossi;i.i, 

• Department of Earth and Planetary Sciences. 54-8 IV. Massachusetts Institute of Technology, ' : 

Cambridge, Mass. 02139. U.S.A. 

( Received 27 April 19721 . 

Abstract — An approximate molecular orbital theory is presented, tested, and applied to silicon and 
aluminum oxyanions. The orbital structure of SiO ,' 4 is calculated and is used to assign the X-ray 
emission spectra. Calculations arc performed at several internuclcar distances and semiquantitative 
agreement is found with experimentally observed trends in AIK<« and Si Kft spectra. Calculations on 
Si s Oj“ ! and AlSiOr r yield fair agreement with experimental trends relating to degree of SiO < poly- \ 

mcri/ation and to Al/Si ratio. Charges for bridging and non-bridging oxygens in Si.O»'\ combined . 
with point charge potentials from metal ions in the orthopyroxene structure yield Ols binding energies :• 
in agreement with ESC A results. 

I. iNTRoni'crioN period elements, ionic models have been used • 

In recent years molecular orbital (MO) to relate experimental line shifts to effective", V 
theory has found increasing application in atomic charges[81, with little use being made 
geochemistry. The need for using MO theory of MO theory since the Ka transition presum- 
in the description of covalent minerals, such ably involves orbitals which are purely atomic, 
as sulfides, has been recognized for some time MO theory has been used qualitatively to as- 
[1J. The application of MO theory to minerals sign Kf3 spectra in such compounds and to 
generally considered to be ionic, such as ox- relate shifts in spectral lines to changes in 
ides and silicates, has only recently become bonding energies. K[3 satellite lines have also; , 
common(2-4]. been interpreted using MO theory. A large; 

A related development is the increased use amount of information on the shift of X-ray' ' 
of MO theory in solid state science generally, lines with average interatomic distance and 
and particularly in X-ray spectroscopy. Sever- with the type of silicate linkage in silicates and 
al authors have suggested that the MO model aluminosilicates is avai!able[9, 10] but has not 
is more appropriate than the atomic or band yet been interpreted theoretically, 
model for the interpretation of X-ray emission - TheonlyquantitativeMOcalculationspre- 
and absorption spectral 5-7]. Qualitative MO viously applied to X-ray spectroscopy are 
diagrams have been used to assign X-ray those of Manne on sulfur and chlorine oxy- 
transitjons and the experimental transition anions] 1 1. 23], The care with which Manne 
energies have then been employed to deter- treated core orbitals and the general ade- 
mine quantitative MO diagrams empirically, quacy of his method resulted in good agree- 
This procedure has been used extensively for ment with experiment for Ka,. t shifts if the 
compounds of the transition elements. relaxation of the core ion states involved in the 

For Ka spectra in compounds of the third X-ray transition was taken into account. Re- 
; cent ah initio self-consistent field (SCF) cal- 

*This work was supported in part by (he National Aero- C Ulations(12| 00 thd sulphate ion have yielded 
nautics and Space Administration. results similar to those of Manne. Although 
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Mannc also 'considered (lie ESCA (Electron from u/> m/7/o SCF calculations on related 
Spectroscopy for Chemical Analysis) or molecules. This method is very similar to one 
photoclcctron spectra (I’ES) and the other presented by Hillier[ 1 5] and applied to SO.," 2 , 
oxyanions data for SO , ' he could not obtain However our method is wholly non-empirical, 
any conclusive results due to a number of un- using non-empirical parameters in calculating 
certainties in both the experimental situation the off-diagonal terms and non-empirical 
and the theoretical model. quantities in the diagonal terms. Also all elec- 

In this paper ,wc present approximate non- irons are included explicitly (and treated 
empirical MO calculations on several, silicate equivalently) and full self-consistency is re- 
compounds. Our goal is to obtain knowledge quired. 

about bonding in silicates and to reproduce To test our method and programming we 
and interpret trends in X-ray and ESCA trans- carried out a calculation on the SO, ~ 2 ion. 
ition energies. We present first the theory. Best atom minimum basis set exponents [16] 
then details of the computation, and finally were used. One center K's were taken from 
o'ur results. atomic SCF calculations and two center K's 

' • for the S-0 bond were assigned values from 

2 . THEORY- an ah initio SCF calculation on PO (phos- 

The present approximate MO calculations phorous monoxidc)[ 17], Calculations on the 
employ a method incorporating elements of hydrides SiH,[l8|. PH : i[17], H 2 S[18] have 
the CNDO/2 method of Poplc and Segal [ 1 3] shown that there is very little variation of 
and the NEMO method of Newton ct <//.[l4]. either one center or two center K fi 's within 
The formulas used for diagonal and off-diagon- such a series. This suggests that a single set 
al elements of the Fock-Hamiltonian (F of F's will givc reasonable results for the ser- 
matrix) arc ' ies Si0 4 PO,‘\ SO," 2 . Also, other SCF 

calculation[19] have indicated relatively little 
F,, = F° + 2 APj {(TiljJ)- . 5(iJliJ)) K variation with charge, electronic configura- 

jonA (1) tion, intcrnuclear distance or orbital expon- 

ent. Results were compared with ab initio 
+ "£ QhR~ a }, minimum basis set [ 1 2a] and various approxi- 

mate calculations [20], The average absolute 
F n = T u + Kj + U jd/2. (2) error of our valence orbital eigenvalues with 

respect to those of [12a] is 0*1 2 AU. very 
where FI is the diagonal element for a free similar to that in Manne's calculation, with 
atom or ion with a specified electronic config- our values being systematically too positive, 
uration. obtained from an atomic SCF calcula- The average relative error in the numerical 
tion. A Pj is the change in the Mulliken orbital ordering of the MO’s is only 0-05 AU. In- 
population of Slater orbital j. (TiljJ) and (ijljj) spcction of the tabulated comparison of vari- 
are rotationally averaged one center Coulomb ous methods given in [20] and reproduced 
an<3 exchange integrals, Q„ is the Mulliken with our addition in Table I shows that our 
charge of atom 8 . (/,, is the rotationally aver- method, as well as those of Manne and Hillier. 
aged element of the potential energy matrix give results similar to the minimum basis set 
and Tji and are the kinetic energy and over- ab initio calculation. Extended basis set cal- 
lap integrals. The first equation is similar to culations yield the opposite ordering for IF, 
that of [13] while (2) is obtained from [I4|. and 5 1,. Correct MO orderings are then ob- 
The overlap, kinetic energy and one center tained except for the lc, and 5 f 2 using all 
electron repulsion integrals are calculated ac- three methods. Although our method yields 
ctirately and the K ti parameters are obtained eigenvalues for ‘core’ orbitals, comparison 
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Table I. Eiitcmalues of SO, - valence mol- 
ecular orbitals using various methods 



Present Ahmiiin Approx, 
method SO-' 1 \tO tw 

sc.c- 

I.CAO 

EHMO'*" 

1/, 

0-371 

: 0-266 

0-152 

— 0-070 

-0-323 

Si.. 

11-434 

0-270 

0-122 

-0-068 

-0-306 

ic 

o :% 

0-204 

0 093 . 

-0 105 

-0-352 

4r. 

0-0 19 

-0-029 

-0-216 

-0-260 

-0 426 

Sit, 

-0-002 

-0-122 

-0-276 

-0-298 

-0-294 

3/.. 

-0-512 

-0 684 

-0 634 

-0 806 

-0 991 

4//, 

-0-718 

-0-933 

-0-848 

— It 18 

— 1-138 


' •'Ref. 1 1 
""Kef: 1 1 1 1. 

'• 'Kef. 1 1 5 1 scll-coiisislcnl chsirjie I CAO. 

‘"Ref. |2»h| extended Htickcl MO. better results have 
been obtained using the nun-self consistent EHMO 
method Ref. I .'I |. 

with the ah initio results is not appropriate 
since they were obtained using Gaussian orbi- 
tals which yield poor energies for inner shells. 
A possible difficulty in our calculations is a 
somewhat exaggerated charge separation as 
illustrated by a S charge of 3-2 in SO,-. 

Computational model and details 
The compounds we discuss in this paper in- 


computation time and severely limit the size 
of systems to be studied. From the ah initio 
calculation |T 2a] and the qualitative work of 
Cruickshank|2 1 ) we would expect the neglect 
of d orbitals to have a strong energetic effect 
only on MO's of c -symmetry, although some d 
orbital participation in the higher /« orbitals is 
indicated by the L spectra of second-row 
oxyanions[22]. 

Orbital structure ofSiOf 4 

Eigenvalues and orbital composition for the 
SiO," 4 unit are given in Tabic 2. The orbital 
structure is similar to that for SO ," 2 . As noted 
earlier our method docs not obtain the proper 
ordering of the 5r. and 1 r , orbitals in SO.," 2 
and this ordering is presumably incorrect in 
SiO,' 1 as well. The inverted ordering of the 
5a, and the 4 1 ., in SiO," 1 is also doubtful al- 
though some shift in this direction may occur. 
The valence orbitals fall into three groups: (1) 
4a, and 3 / 2 : 02, (2) 4t 2 and 5a,: a bonding 
orbitals: (3) \e, 1/,, 5/ 2 : non-bonding orbitals. 
The electronic configuration expressed as 
Mulliken orbital populations is Si(core)' n 
3.v o;w 3 p' ,m O(core) 2 2r 2n:, 2p SS2 silicon and 



volve Si coordinated with oxygen. We first oxygen. 

present calculations on the perfect tetrahedral * ’fi'* 

anion SiO r' for various Si-O distances. A X-ray transitions ^ 

best atom minimum basis set (BA MBSf of The quantitative MO results of Table 2 J 

Slater orbitals is used and K s arc again taken ma y he used to assign transitions in the X-ray 

from atomic SCF results and from PO. Cal- 


culations were assumed to have converged Table 2. Eigenvalues and orbital composi- 
when the average F„ variation was < 0 005 ,ions of Jilted orbitals for SiO f 4 (R(Si-O) = 
A U. This required 4 or 5 cycles. 1-603 /4) 


Some controversy has surrounded the use 
of 3 d orbitals in MO calculation* on third 
period compounds. We have chosen to neglect 
3d orbitals for the following reasons: ( I ) com- 
parison with accurate SO,"'- calculations 
showed our method to give significantly poor- 
er agreement when 3d orbitals were included. 
(2) use of a 3d function with a minimum v. p 
set exaggerates d orbital participation. (3) d 
orbitals are not needed for an adequate des- 
cription of most properties of SO,'-’. (4) an 
increase in the basis set size would increase 


i 

edAU) 

■ " % composition 

St, 

1 -07332 

Iff- Si3/7. 99% 02/7 

1 

1-02677 

t00%O2p 

u- 

0-96229 

100% 02/7 

5«, 

0-85682 

!7%Si3.r. 1% 02*82% 02/7 

4 1 . 

0-84077 

19% Si3/7. 7 \% 02/7 

3 f. 

0-18179 

3%Si3/>. 97% 02.1 

4«, 

0 11357 

3% Si3r. 97% 02.v 

2 1. 

-3-22745 

99% Si2/7. 1% Si3/7 

•Vi 

-5-04019 

1005? Si2.v 

1/. 

-18-67688 

100% Ol . v 


-18-67705 

I00 r ? Ol.v 

!</, 

-67-98569 

100% Si 1 v 
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rfpcs r dn Pif. 

spectra of Si and O in silicates. MO theory SiK B 

has previously been used for this purpose by The Si At 0 spectra nominally involve the 
Urch |6J and by Andermann and Whitehead transition "ip -* l.r. In fact several bands in 
|2.H X-ray spectra are produced by transi- the K„ region arc observed and by the dipole 
tions from a higher level to an inner shell va- selection rules for tetrahedral geometry all 
cancy and may be classified as AC. L or M molecular orbitals of t, symmetry can take 
depending upon the inner shell involved. The part in the transitions. In Si0 2 three such 
AC„ line arises from the transition 2/> — l.v and transitions arc observed, the main transition 
the SiACa energy is thus given by the difference being called K tt and the lower energy transi- 
of the 2 1, and Id, eigenvalues in Table 2. Due lion K' ir Table 3 contains experimental and 
to the approximate nature of the calculation calculated energies and intensities relative to 
and to the neglect of the important relaxation AC„ for Si0 2 . The AC 0 and K' & experimental val- 
effccts associated with core ion states agree- ucs arc taken from Urch [6] while the position 
ment of the calculated absolute value of transi- and approximate intensity of the higher energy 
tion energy with experiment is only approxi- peak arc obtained from the figures of Dodd 
mate. For example; we calculate £' S nn as and Glenn. Relative intensities are calculated 
~ 1761 cV while experiment gives ~ 1740 cV. in a simple way from the relative 3 p pcrcent- 
• (see .the Appendix for further discussion), age compositions of 3/ 2 , 4/ 2 and 5/ 2 in Table 2, 
Therefore we will consider only the relative thus ignoring ‘cross-over’ contributions (e.g. 
energies of the spectral transitions. 02 p -* Sils) which Urch has shown to give a 


Table 3. Experimental and ealeuhitcd energies < cV) and intensities 
of. SfK*. vX/ls.*- and 'OK„ spectral components for tetrahedral 

Si-O coordination 



Expcr."" 


SiAC„ 

Calculated'" 1 

k;, 

±1- 

t 

Assignment 

AE 

1 ' 


— 13-7. — 14-4 

15-20 

3r..— In, 

-17-9 . 

16 

A n 

0 

fiio 

4 1 , — In, 

0 

100 

shoulder 

+5. . 

- 1-5 

5/- — In, 

6-3 

• 5 


Expcr. 1 ' 1 



Calculated 


At: 

/ 

Assignment 

AE 

/ 


0 

20-25 

4d, — It,. 

o 

18 


• i: 

100 

5u, — 2/.. 

20-2 

100 


- 17 

90 

Ic. 5 — 2/.., .. 
OK a 

22-9, 25-7 

0,42 


Exper."" 



Calculated 


AE 

/ 

Assignment 

AE 

/ 


0 1 


5/.. — (l/-.2<»,) 

0 

99 


-10 

-:-8 ■ 1 

-4 

1/, — * ( 1/-.. 2(/,l 

-1-2 - 

100 



• Ic — (U..2u,> 

-30 

67 


-y* 

1 

5((, — (1/-. 2(/,l 

-5-9 

27 




4( . ( 1 ( 2((,l 

-6-3 

71 


'■''SiO, U 1-603. 
"•'Ref. |2|and |f*|. 
'"Rcf.|22|. 
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small intensity. Manne|24| has shown this 
approach to give good results for hydrocar- 
bons. Agreement of calculation and experi- 
ment is fairly good. We do not find it neces- 
sary to invoke a symmetry forbidden transi- 
tion from the 2c as did Dodd and Cilcnn|21 in 
order to explain the high energy shoulder. Our 
assignment is supported by the results of cal- 
culations using the same method and BA M BS 
exponents on isoclectronic AIO,~\ The cal- 
culated separation of main K u peak and high 
energy shoulder is about i -6 eV less in AI0 4 _i 
than in SiO,~\ in agreement with the differ- 
ence of 1 -5 eV in separation observed be- 
tween quartz and cordierite[32). 

OK„ 

OK a spectra arise from 02p— >01S 
transitions. The upper orbitals contributing 
to OK a are 4t... 5a,; I c. I/, and 5t->. Table 3 
gives calculated energies and intensities for 
the O K„ spectrum obtained in the same way 
as for the Si/f u . In this table thc highest energy 
component is assumed to arise from 5 t- : al- 
though the It, may well lie higher. Agreement 
with the K n structure as discussed by Urch is 
fairly good. A broad peak is found resulting 
from transitions from the non-bonding pre- 
dominantly 02p orbitals and a lower energy 
peak results from the <r bonding orbitals. The 
assignment of Urch also mentions an expected 
line at about 4 eV lower energy which has not 
been observed. Urch assigns this observed 
peak to the 5a,. but we find 4t._, and 5n, to be so 
close in energy that they would generate only 
a single peak. The calculated intensity ratio 
for the two bands in the K„ spectra, summed 
over ‘components, is approximately 2-7/1. 
somewhat smaller than experiment. We note 
that our MO diagram is not in agreement with 
that of O'Nions and Smith 1 3 ). We feel that 
the high energy peak at 531 eV in |3| cannot 
be accomodated on a conventional MO dia- 
gram; it probably arises from doubly ionized 
aloms|24|. On the other hand if this peak is 
ignored, -agreement of our calculations with 


O N ions and Smith's reported spectra is good, 
as shown in T able 3. 

SiU, 

SiL-. 3 . spectra nominally involve the transi- 
tions Si 3x. 3d -* Si2p. The upper levels in- 
volved arc 4«,. 5fl| in our MO scheme, both 
having some Si3s character. Clearly d orbital 
participation is required in this case to give 
the three or four strong lines observed. 

Using our SiO., -4 results in conjunction with 
the ab initio s, p. d results of Gelius et al. for 
SO.,' 2 We can suggest the assignment of the 
peaks which is given in Table 3. We assign the 
low energy peak to 4a, and the middle peak to 
5a,. The peak at higher energy must then be 
assigned to lt» and/or 5 t t . In the ah initio cal- 
culation on SO,~ 2 substantial d orbital popula- 
tions were found in le and 5 1 ,. Our interpreta- 
tion is also generally consistent with an ah 
initio SCF calculation on Si0 4 -4 , currently in 
prcss(25]. 

Earlier MO assignments of the L spectra of 
S and Cl oxyanions ignored intensity con- 
siderations and claimed agreement of calcul- 
ated L spectra with experiment using an v. p 1 
MO calculalion[23]. However, although d 
orbital participation may not drastically affect 
the eigenvalues of the \e and 5c orbitals it 
will drastically aficct the intensity of the L 
spectral lines generated by them and so an s.p 
basis cannot yield both correct energies and 
correct intensities for Si L 2 . ; , spectra. 

This is seen in the L s<3 intensities calculated 
in Table 3. The intensity ratios of the two 
lower energy transitions were obtained from 
the %Si3.v character in 44, and 5a,. To obtain 
any intensity for the high energy peak we had 
to consider the elements of the dipole mom- 
ent operator connecting the Si2p and 02 p 
prbitals. These integrals were calculated ac- 
curately and we found the 02 per electrons of 
5 1 2 to generate a transition ~ 42 per cent as 
intense as that from the 5a,. The intensity 
from the 02 ptt orbitals of l<* was essentially 
zero. Our L.., intensities are then in poor 
agreement with experiment, although the 
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energies are. quite good. The poor intensities 
are clearly due to the neglect o f </ orbitals in 
our calculation. 

Experimental amt calculated variation oJ'E sili 
and Estk with R(Si-O) 

Calculations at the level of the present work 
arc more appropriately applied .to the study 
of trends in a scries of similar compounds than 
to the calculation of absolute values of a 
quantity in a single compound. Therefore wc 
have studied the variation of K spectral ener- 
gies with variation of several different experi- 
mental quantities. Wc first varied R(Si-O) 
for tetrahedral SiO , -4 . Trends in E hts arc dis- 
played in Fig. Kb). The variation in E hn is 
about 4-5-5 x 10";' AU/0 01 A. Studies of S\K n 
spectra by Koffman and Moll 1 26] have shown 
no disccrnable variation in £' SUa over a large 
range of compounds. On the other hand 
Wardle and Brindley [10] have found fc' A uv„ to 
increase as /?(AI-0) increases in both tetra- 
hedral and octahedral coordination. Their 
experimental slope of L vs R{ Al-O) for A1 
in Td coordination ((KAl-O) — 1-7-1 -8) is 
3-4 x I0" 4 AU/0 01 A. We performed calcula- 
tions for Aid,"', analogous to those for iso- 
electronic SiOi" 1 . at ft(AI-O) = I -70. 1-75. 

1 -80 and obtained a slope of 5 x 10* 1 AU/0 01 
A so that our agreement with experiment is 
good for E a | A „. Experimental and calculated 
Al/f« results are given in Fig: Ha): The Si AT„ 
results will be discussed further later in this 
paper. 

At this point we note for convenience that 
in the Si K 0 energy range a variation of I X 10“' 
A in wavelength is closely equivalent to an 
energy difference of 0-001 AU. In the work of 
White and Gibbs on K„ shifts in silicates | Va| 
5iAf„ shifts were given relative to SiO_. (or- 
quartz) in units of I x 10 ‘ A (their A scale). 
They found E s iK to increase rapidly as W ( Si— 
O) av increased and fitted a straight line to 
their plot of A vs /KSi-O). Wc feel that their 
plot of A vs W(Si-O) fails to distinguish be- 
tween two separable, although related, fac- 
tors which affect A;' Si/k ( I ) ,/HSi-O) and (2) 






t iii. I. la) l-xpcrimcntal data from (10| compared 

with, calculation — experimental calculation. 

(h) I xpcrimcnlal data from |‘)a| compared with cal- 
culation. t he zircon value of ft(Si-O) is 1-6222 A. the 

recently relincd value of 1 301. new isolated tetra- 

hedra line calculation original line from |9u|. 
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the degree of polymerization of the silicate 
tetrahedra. We have therefore rcplotted their 
data for their isolated tetrahedra minerals only 
and obtained a line uith a slope of — 1-4 A / 
0 01 A (0 0014 AU/0-01 A), much smaller 
than the slope of their original line. Over the 
range R{ Si-O) = I 043-1 003 A we calculate a 
slope for K„ of 1 • 1 A/0-01 A. in good agree- 
ment with the new line. This is shown in Fig. 
Kb). 

Since Al ,:1 minerals arc substantially more 
complex than those of Si '. at the present time 
we must generally ignore them. However, we 
note that the slope of I: M k„ with /KAl-O) in 
octahedral coordination is the same as that of 
F.Si K„ with R( Si-O) for silicon in tetrahedral 
coordination|9b|. We would expect Al or Si 
in octahedral or tetrahedral coordination to 
behave similarly and we have found similar 
K„ and K„ trends in calculations on Si and Al 
oxides in tetrahedral and octahedral coordina- 
tion. 

(dlenlatcd variation of ( with decree of 
polymerization and with A II Si ratio 

An accurate determination of the 
shift with degree of polymerization is impos- 
sible at this time, since very large polymeric 
systems would necessarily have to be con- 
sidered. We can however perform a calcula- 
tion for the smallest fragment of the linearly 
polymerized ( SiO ; ,)„ chain, i.e. Si.CV'f £W* 
values are given for several geometries in 
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Table 4. For these calculations we generally 
find./fsu to be lowered from its SiO^" 4 value. 
In this C‘2, system seven non-degenerate '4t- 2 - 
likc* orbitals contribute to the transition and 
we must weight their eigenvalues by their per 
cent Si3p character to obtain an average £ siAo . 
1 he energy splitting of these ‘4/ 2 -like’ orbitals 
in Si..O ; _,i is substantial. The energies of the 
resulting K n spectral transitions are shown in 
Fig. 2 along with the %Si3/> character of the 
orbitals. A later paper will discuss the splitting 
of Ihc ‘4/^-like’ orbitals in more detail. The 
variation of £ si/l/J with 0 = L Si-O t)r -Si is 
quite simple. Figures 3(a) and (b) and Table 4 



Table 4. Values of!-j Sth and related quantities for SiO t ~\ SLOf 11 and AlSiO-f 


SiO, 

• 

K ■ 16(0 

• i 

K>U 

Si,<> : 

i) ■ 120 

141- 

K(Si-O) 

160° 

= 1 603 
180’ 

Si, Or" 

R(Si-0„ hr ) == 1-603 
R(Si-0„,)= 1-623 

AISiOr‘ r 

4 '• 

0 004 

110(166 

-0-0033 

-0 0057 

-0-0018 

-0-0030 

0-0366 

It . 0 

00022 

(10033 

0-0035 

' 0 0021 

0-0010 

0-0026 

0-0158 

DiSii 2-582 

2-328. 

2-666 

2 664 

2-664 

2 653 

2-631 

2-664 

m Si •()..,) 

— 

ouoi 

0-096 

0-097 

0-090 

0-096 

0-196 

1/2 v „ jf/ 

— 

- 268 64 

-269 (IS 

-269 14 

-268-88 

-269-19 



Pncnjic'. .ire in A l . 

m St— < 1 1 v ;i imillikcn overlap population between Si and <) hr . 
\ll . , -Si-O,,,., lire 





12 


ft sea, m w 


J. A. TOSSlil.l. 


-C 



0.08 




0.06 


t; 


„. 0.04 

u 



° 0.02 
<5. 

b 

v 


> 

^ 0 
CO 

« -0.02 

i 

J 


-0.04 


u 

C 

-0.06 

-0.08 


< J , 
h- ^ 


$ 



Kn 


120 140 160 180 

8 
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Fig. 3(b). Variation of V , ( ,e, with tf (AU). 


show that £ S | A<j : varies with 0 in the same way 
as does the sum of the MO eigenvalues. Both 
quantities have a minimum near 0 — 160° when 
all Si-0 bond lengths are constrained to be 
1-603 A. For 0 = 160° the lowering of £ S iA fl is 
about lending support to the hypothesis that 
the degree of polymerization is an important 
variable influencing £ MA/j A very approximate 
formula for the K u shift of a silicate relative 
to «- SiO- can then be given (in A units) as 


A/: S1A(J = 100 { R{ Si-O)^" 1 " 1 ' - I -607M 


+ (4 -/»“')« < 3 > 


where li*"™"-' is the number of oxygens in- 
volved in bridging bonds in the silicate. Our 
calculations yield A and B values of about I 
and 3-6 respectively. The correlation of the 
eigenvalue sum with the eigenvalue of the 4f 2 
valence orbital is expected and is consistent 
with the interpretation of Dodd and Glenn [2]. 
To determine the relation of the eigenvalue 
sum to the total energy would require further 
calculations since the total SCF energy is 
given by the expression 

E = E sn + 2 +'Zn,H i (4) 

i I 

where E y „ is the nuclear-nuclear repulsion 
energy and Hi is the one electron energy (// = 
T+V) for MO i. Assuming that eigenvalue 
sums parallel total energies for variations of 8 
the results of Table 4 are consistent with the 
known relationship between /f(Si-O b r) and 8 
[4]. Bridging oxygen bond lengths are general- 
ly larger than those for non-bridging oxygens. 
As fl(Si-Ohr) decreases toward non-bridging 
values 0 increases. Therefore, if all Si-0 bond 
lengths are set equal we would expect 8 to 
increase from the value of 141° which results 
from linkage of two SiOr 4 tetrahedra without 
any O-Si-O angular distortion. The calcul- 
ated minimum energy value of 160° is quite 
reasonable. In later calculations we hope to 
explore the relationship between bond length 
and 8 in more detail. 

We find E S \K a to increase by about 0-003 
AU in going from SiO, -4 to S^Ot" 6 , consistent 
with the higher Si charge, although £ S u 0 de- 
creased as the silicon charge increased for R 
variation in SiO, -4 . This suggests an explana- 
tion for the different behavior of £* for Si and 
Al. For Si the variation of /?(M-0) is quite 
small and smaller values of R are associated 
with increased polymerization. The effects of 
reduced 7? (Si-O) and increased polymeriza- 
tion on £' S ia are in opposite directions and 
<* 

tend to cancel. For Al the /?(M-0) variation 
is considerably larger and its linkage to poly- 
merization is more complex but generally in 
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the opposite direction, so that a clear trend in 
I£ .\\ k results. For F sjA(j the trends are additive 
and therefore substantial i£ K variation is ob- 
served. fo 

A final trend observed in Si/v'„ spectra is 
an increase in E SIAfl with an increase in Al/Si 
ratio. For a AISi0 7 ~ 7 unit with ft(AI-0) = 

1 -70 and /?(Si-0) = 1-603 we find £ SIA(j to be 
raised in energy by about 37A with respect to 
SiO,' 4 . The variation found by White and 
Gibbs between framework aluminosilicates 
with Al/Si = 0 and those with Al/Si = I is 
about 25A. Gur calculations therefore some- 
what exaggerate the effect of Al substitution 
but give a trend of reasonable magnitude. The 
mixing of the Si3/> orbitals with the higher en- 
ergy AI3 p orbitals is the reason for this result. 

Some experimental data on Si K,, spectra 
also exists for a high pressure polymorph of 
SiO.- stishovite, with Si in six-fold coordina- 
tion in a framework structure. Calculations for 
the octahedral species SiOr,"* will be discussed 
in detail in a later publication.^ We note here 
that C.siA w for, SiO„ _H at /f(Si-O) = T80 A is 
higher than for SiO.,' 4 by about 70 A. This is 
in fair agreement with the experimental 7;*^, 
difference of about 50 A between n-quartz and 
.'stishovite and suggests that differences in K a 
spectra arising from differences in coordina- 
tion number may be reproduced fairly well by 
our method. 
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Charge distribution information from SiK a:Ui 
shifts 

The SiK, tl A satellite lines from doubly- 
ionized atoms can also yield information on 
electronic structure. Although Koffman and 
Moll (26) stated that there was no obvious cor- 
relation between their measured S\K a energies 
and chemical nature in the group of silic- 
ates they studied, we find that their K „ 3 , 4 
data fall naturally into four groups: Si0 2 , 
other framework silicates, chain silicates, 
and isolated tetrahedral silicates. The 
results are consistent with Si charge which 
increases with the degree of polymeriza- 
tion. In Table 5 we present our grouping 
of the K„,i. i data of Koffman and Moll. From 
isolated tetrahedron ( Fe,SiO,) to frame- 
work structure (SiO..) we find a decrease in 
K hn: of about 7-8 x 10' 4 A. Using the value 
for the variation of with atomic charge 

given by Chun [8], we calculate a difference in 
effective Si charge between the two silicates of 
about 0-25 electrons, with Si more positive in 
SiO,.. The calculated <2, sn values in Table 4 
show a trend in agreement with these results 
both in direction and magnitude. 

The results of Dodd and Glenn give very 
high values of A for some pyrosilicates 
(SioOj ” groups), inconsistent with our results. 
However these minerals are unusual, e g. in 
one of them, gehlenite. there are three T d posi- 
tions with R(M-0)= 163, 1 -63 and T87 A. 


Table 5. Regrouping of the \SiK a .,., (A) data of Koffman and 

Mol!\ 26J 


Group 

K a., 

Ktt-% group 
average 

Ka, 

Ka, group 
average 

Group 

structure type 

SiO. 

7-0765 

7-0763 

7-0673 

7-0673 

Quarts 

KAl'Si.O. 

7-0765 


7-0678 


Framework 

NaAISi.O, 

7-0766 

7-0765 

7-0677 

7-0677 

aluminosilicates 

CaAl.Si.O. 

7 0763 


7-0677 



Mg.Si-6, 

7-0768 


7-0680 



• CaSiO, 

7-0770 

7-0770 

7-0680 

7-0680 

Chain silicates 

KAISi/),, 

7-0772 


70680 



CaMgSijO,, 

7-0770 


7-0682 



Fc,SiO, 

7-0773 

7-0773 

7 0680 

7-0680 

Isolated tetrahedron 

Si metal 

7-0803 

7-0803 

7-0713 

7-0713 

Metal 


increases as degree of linkage inercases. 
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Thc distribution of silicon between the three 
sites is not known and so the e licet of/ftSi-O) 
variation cannot be determined. This may ex- 
plain the discrepancy. z 


following formula: 

E^ =. i + s A /y{ (( /7/U) - 0-5 Okiik)] - 

[(Jjlkk) — 0-5{jkijkj ] } (5) 
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Predicted shift of other lilies 

Our results allow us to predict the shifts in 
other X-ray emission lines. The Si K' B line 
should shift to higher energy with increasing 
R at about three times the rate of K. However 
with increasing polymerization /: SjA/) . increases 
so the effects tend to cancel and we expect 
little /:' SiA . variation. In the Si L..-* spectra wc 
associate the peaks with the 4o,. 5 n, and 5f 2 
orbitals. As R decreases and polymerization 
increases, the effects of the variations again 
cancels so these peaks should show little shift. 
For the OK a peaks there is a general decrease 
in energy as R decreases, but A is higher 
for O hr than for O n ,, r in Si.-O? - ' 1 . We feel that 
the effect of polymerization will predominate, 
leading to higher link as the degree of poly- 
merization increases. Although such experi- 
mental comparisons arc not available for sili- 
cates we can compare E l)h for SiO , and Al ,0 : , 

Ct 

|27]. The OK a line is found at lower energies 
in Af.O.i by about 0 08 AU. This is consistent 
with out calculations which show O to be 
more negative in AIO.» -3 than in SiO, -1 (the 
difference is about 0 07 electrons). For the in- 
equivalent oxygens of SbO?;" we find E 0A to 

Of 

shift to lower energy as the oxygen becomes 
more negative at a rate of T5 AU/electron. 
so the calculated difference of O charges in 
SiO. and AUO ;! could yield the experimental 
result. 

Simplified model for SiK shifts 

Employing equations (I) and (2) as a basis, 
wc can develop a, simple model to explain SiK 
shifts. As noted by Mannc the energy separa- 
tions of core orbitals can be obtained from one 
center integrals alone, since the two center 
terms cancel (even in accurate calculations). 
We represent the dependence of the energy 
separation on the orbital populations by the 


where is the transition energy in a refer- 
ence state and A P k is the difference in k orbital 
population. We can abbreviate this as 

and in Table 6 wc give values for S,. a - B . We see 
from Table 6 that while the effect of 3.r and 3 p 
populations on E hn will be appreciable the 3 d 
population will have little influence. For this 
reason calculations with and without 3d orbi- 
tals. giving different Si charges but similar 3.r 
and 3 p populations, would probably give simi- 
lar K„ shifts. The term involving 8"- ( has norm- 
ally been ignored, but this is an important term 
influencing the variation of £ Aq very strongly. 
It should be remembered that this 2p orbital 
is a Slater orbital rather than an atomic Har- 
tree-Fock 2 p orbital. Equation (6) has been 
found to give good results for E K for all the 

a ; « 

systems considered in this paper, verifying 
the opinion that bonding effects (as opposed to 
charge distribution effects) may be neglected 
in determining £si/v u - 

The values of 5° are much larger than those 
of 5". consistent with the larger dependence of 
£ A/j on charge found in atomic SCF calcula- 
tions. However since the Si3p is a valence 
orbital equation (6) is not accurate. A satis- 
factory simple formula-for £ Si * fl has not been 


Table 6. Value's of 8*° 
andb K ° (AU) 


K 

V 

6/ 

3i ,#l 

0 0020 

0 1525 

3 p u " 

00138 

02917 

yd "" 

0 0006 

00826 


1 5513 

1 9922 


'“•’B A MBS exponents. 
"•'|=l-3;, , , 
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found. The main reason seems to be that the 
‘4/,,-like’ orbitals have both Si 3 /j and 02/) 
character and the dependence of their eigen- 
values upon the /-' matrix elements is quite 
complex. In genera! increased SB/). 02p 
overlap is associated with lower Si AT 0 values 
for the SiO.j" 4 series. For 0 variation in Si_.0 T ~ 6 
£ s1A(( correlates with the eigenvalue sum 
which should roughly parallel the total mole- 
cular energy. 

Discrepancies between calculation and experi- 
ment for Si X-ray emission spectra in Si0 2 
and Si metal 

Experimental data also exist for the Si K„ 
and K a spectra of Si in silicon metal. In com- 
paring SiOo with Si metal, experiment shows 
£ sl * n to increase by 0 016 AU in going from 
metal to oxide and E S \k ii to decrease by 0160 
AU. We find E S1A - for SUV' 1 to lie at 0 072 

.a • 

AU above the free Si atom (2p- IS) value. Al- 
though the Mulliken charge of Si must be zero 
in Si metal, we would expect the valence elec- 
tron density to become mbre diffuse and E ha 
to increase as though the Si possessed a partial 
positive charge. This explains part of the dis- 
crepancy, but another part must be due to an 
exaggeration of the Si-O charge separation in 
our calculations. 

We find the E st K g for ShO? - " falls far below 
the line determined by the one center term of 
equation (6) alone, but is 0-27 AU above the 
free atom value since £V fl increases rapidly 
with Si charge. The explanation for the dis- 
crepancy must be similar to that in the K a 
case above. The K a and K u energies of Si TV 
would both have the proper "Values with 
respect to Si metal if the difference in effective 
Si cltarge between Si-jOr'* and the metal were 
about I I. The large lowering of E Kg from the 
value given by equation (6) is essential to this 
result. This lowering is clearly associated with 
Si3„. 02,, overlap. Studies of transition metal 
K„ : , spectra have lead some. authors to calcu- 
late effective charges' of such atoms in the 
metal from the variation of with oxida- 


tion state 1 5]. Charges around + 1 are generally 
found. If Si is assigned an effective charge of 
about + ! in silicon metal, our results are in 
reasonable agreement with experiment.. • 0 

s ko~lc/ be 

ESCA shifts Ty^’k-r’o. . i 

The shift of the Ols ESCA line in silicates 
has been studied by Huntress and Wilson (28] 
and by Yin, Ghose and Adler[29]. Huntress, 
and Wilson found the Ols binding energy 
(BE) to be lower in fayalite (Fc ; Si0 4 . isolated 
tetrahedra. all oxygens non-bridging) than in 
SiO> (framework, all oxygens bridging) by 
0-5 eV. Yin er «/., on the other hand, find non- 
bridging oxygens to have higher binding ener- 
gies than bridging oxygens in orthopyroxenes. 

Such apparent inconsistency demonstrates 
that c' taining information on charge distribu- 
tions from solid state ESCA results will often 
be difficult. In our Si.,0 7 ~ 6 results we find O b > 
to be: less negative than O nbr by about 0- 13 
electrons. By Koopmans theorem [34] we 
equate ionization potentials and eigenvalues. 

The e (Ols) difference between Q br and O nbr 
is about 0-21 AU with O br Is more strongly 
bound. To calcalate relative Ols BEY in 
orthopyroxene we must aiso consider the ef- 
fect of the metal ions needed to produce 
charge balance: It is well known that the effect 
of these ions may be large[ll], [35). In py- 
roxene the distances of the bridging and non- 
bridging oxygens from the metal atoms are 
very different. The metal in the more distorted 
site (designated M2) is in 6-fold coordination 
with 4 O nbr arid 2 O br with the distances to O bt 
being almost 0-5 A greater. The metal in the 
more regular Ml site has only O nbf in its first 
coordination sphere. Distances from Ml to the 
next nearest neighbor O hr are at least T3 A 
greater. Using the X-ray diffraction data of 
Ghose[33] we have calculated the point 
, charge potentials produced by two +2 charges, 
one in M I and one in M2, at the positions of 
O br and O nbr . Point charge potentials from all 
more distant atoms were ignored since they 
should approximately cancel for O br compared 
to 0„ br . We find the Ml and M2 point charges 
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to stabilize O n) , r by about 0-26 A U with re- 
spect to 0| >r - The point charge effect of the 
metal atoms is them approximately equal and 
opposite to the eigenvalue difference in free 
. SkO;" 6 . Using free Si.0 7 ~ ,: with the point 
charge term added we find 0„ r Is less strongly 
bound than O nt , r Is by about I eV. in fortuit- 
ously good agreement with the results of 
Yin el til. considering the numerous approxi- 
mations in our calculation. 


CONCLUSIONS 

Our results indicate that molecular orbital 
calculations at the present level can be useful 
for assigning X-ray emission spectra. They 
can also yield a great deal of semi-quantitative 
information on the relationship of SiK a and K u 
transition energies to the molecular and elec- 
tronic structure of silicates. Calculated charge 
distributions arc also found to be in agreement 
with FSCA results. Comparison with X-ray 
emission and ESCA data appears to be an 
excellent way of testing the accuracy of MO 
calculations. The reasonable agreement of our 
calculations with experiment suggests that 
chemical trends arc given correctly, particu- 
larly trends in charge distributions. Our calcu- 
lations indicate clearly that Si becomes more 
positive as /? (Si— O) decreases and as poly- 
merization increases and that bridging oxy- 
gens arc less negative than non-bridging 
oxygens. 


Table 7. Experimental and calcu- 
lated X-ray transition . energies 
(EV) 



Exp. 

A 

Calc. 

ft 

C 

K„ 

1740 

1761 

1757 

1732 

...A',, 

1832 

1872 

1858 

1837 

'II. Ill 

89 

III 

1(H) 

98 

OA„ 

526 

531 

539 

524 


.-I uncomxtcil. 

It correct inn ti Knitted 
C enrieelitinc (it .mil (ii Knitted. 


APPENDIX 

Calculated absolute values of X-ray transition energies 
arc in fair agreement with experiment if corrections are 
made for ti) the difference of minimum basis set and near 
I lartrce-Tock basis set atomic SOF eigenvalues and (ii) 
the electronic ‘relaxation' of the core ion. In Table 7 we 
give experimental energies, and our calculated energies, 
uncorrccted. and with corrections of type (i) and (ii) 
added. The eigenvalues of Si*- and O' from our minimum 
basis set and the near Hartrce-Foek basis of Clemenii 
( x ft | form the foundation for correction (i). The approxi- 
mate formulas of Snydcr|37] were used for correction (ii). 
Both upper and lower levels involved in the transition 
were corrected for relaxation. 
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Interpretation of K X-ray emission spectra and chemical 
bonding in oxides of Mg, A1 and Si using quantitative 
molecular orbital theory 
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Abstract— -Quantitative molecular orbital calculations are reported for Mg, Al and Si in tetra- 
hedral and octahedral coordination with oxygen. These calculations are employed to assign and 
interpret the MK a , MK, and OK* X-ray emission spectra of the corresponding oxides. The 
interpretation of tho MR. spectrum reproduces the observed trends in main peak and satellite 
energies with variation of metal, ligand and -geometry. The splitting of the main K p peak, ob- 
served in many oxides, is found to bo a result of interaction between adjacent metal atoms. The 
calculations also reproduce the observed trends in OK a spectra. The electric structures of the 
various oxides arc discussed briefly. 

Introduction 

Molecular orbital (MO) theory is presently finding increased application in geo- 
chemistry. For example, several MO studies of bonding in silicates have recently 
appeared (Gibbs et al., 1972; Louisnathan and Gibbs, 1972). Qualitative MO 
diagrams have been used for some time to assign X-ray emission spectra (Dodd and 
Glenn, 1968). Also, experimental X-ray transition energies have been used in con- 
junction with spectral assignments from qualitative MO theory to construoUquanti- 
tative MO diagrams empirically (Fischer, 1971). For the oxyanions of sulfur and 
chlorine approximate self-consistent-field linear combination of atomic <orbitals 
molecular orbital (SC'F LCAOMO) calculations are available (Manne, 1967) and have' 
been used in the interpretation of the X-ray emission spectra of these compounds 
(Andermann and Whitehead, 1971). 

Since X-ray emission spectra yield well understood experimental quantities which 
.can be easily compared with the results of MO calculations they provide an excellent 
method for testing the validity of such calculations. The different series of X-ray 
emission spectra (e g. SiK,, SiK„, SiL 2 3 , OK, in Si0 2 ) can be combined to give in- 
formation on the energies and charge distributions of almost all of the molecular 
orbitals of the sytem. These experimental results then allow us to assess the accuracy 
of the calculated quantum mechanical results for the whole set of MO’s rather than 
for just a small subset (as is the case for most spectroscopies). Once agreement with 
the experimental X-ray emission data has been demonstrated we can have some 
confidence in our ability to describe accurately, using MO theory, those properties 
which are more difficult to determine experimentally. Examples of properties of this 
type are the binding energies of various geometrical configurations of mineral 
structures and the charge distribution in these structures. In this paper, MO results 
for oxides of Si, Al and Mg are presented and these results are used to assign and 
interpret data from X-ray emission spectroscopy. After agreement with experiment 
has been demonstrated a brief discussion is given of the calculated molecular charge 
distributions. 
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Theory 

The molecular orbitals arc obtained as solutions to the quantum mechanical 
Schrotlingcr equation in its Hartree-Fock-Roothaan form, suitable for SCF LCAO 
MO calculations (Roothaan, 1951). The equation to be solved is, in matrix form; 

FC = eSC (1) 

where F is the Fock Hamiltonian matrix; C is the vector of coefficients of atomic 
orbitals (AO’s) in the MO’s: S is the overlap matrix; and e is a diagonal matrix 
giving the MO eigenvalues, which correspond to vertical ionization potentials. The 
F matrix, the essential quantity in equation (1), is rather difficult to calculate. 
Several schemes have been developed for approximating its elements. The method 
employed here is a refinement of the non-empirical molecular orbital (NEMO) 
method of Newton et al. (1900). It is a self-consistent -field version of this method 
which incorporates elements of the complete neglect of differential overlap method, 
version 2 (CNDO/2) of Pople and Skoal (1900) and is called SCF NEMO. 

In both the NEMO and the SCF NEMO methods the off diagonal elements of F 
are given as 

'Oii + Ojj\ 


E,, = 


T it + KijSij(— 


(2) 


where Tij is a kinetic energy integral, Oii is an averaged potential energy integral 
and Kij is a parameter chosen to reproduce Fij values obtained from accurate ab 
initio solutions to the Hartree-Fock-Roothaan equations for related molecules, in 
which all integrals are evaluated accurately and no approximations (beyond those 
leading to the Hartree-Fock-Roothaan equation) are made. Studies of the theoretical 
foundation of Extended Huckel MO (EHMO) theory by Newton tt al. (1960) have 
indicated equation (2) to be superior to its Extended Huckel analog 


Fij - KSij{™±££) 


For the diagonal elements of F the NEMO theory gives 

Fii 


at 


( 3 ) 


( 4 ) 


where the at are the diagonal elements from ab initio calculations on related molecules. 
In the SCF NEMO theory the diagonal elements depend upon the charge distributions 
in the molecule studied according to the formula (for a closed shell electronic state) 

(5) 


Fii + Fii + 2 A Pj{{iiljj) — 0-5 (ijlij)} 

ton A 


+ 2 QbRab 1 

BsM 


where Fii is a reference value obtained from an atomic SCF calculation for a specified 
electronic configuration of the atom or ion, A Pj is the change in the Mulliken popu- 
lation* of atomic orbital j, Q n is the Mulliken charge on atom B and ( iiljj ) and (ijlij) 

* Definitions of the various charge distributions quantities are given in Mulliken (1955). 
The defining equations are (for closed shells): population of orbital j: Pj — 2 2 2 CnkPkj 

• n k 

charge on atom B : Q B — Z D — 2 P k 


overlap population of A and B 

n(A - B) 


k on B 

2 2 I I 

«i ion A tonB 
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are averaged values of the integrals of the electron repulsion operator e'-r,, 1 , where 
r„ is the distance between electrons 1 and 2. 

Earlier comparison of SCF NEMO with ah initio SCF results for the S0 4 2 anion 
has shown SCF NEMO to give relative MO energies to an accuracy of 1-2 eV,. 
(Tossell, 1972). The SCF NEMO method includes all electrons so that energies are 
obtained for both valence and core MO’s. 

MO Structure of Si Oxides in Tetrahedral (Td) and Octahedral (O h ) 

Geometries 


The gross MO structure is found to be the same for both the tetrahedral and 
octahedral oxides of Si, A1 and Mg. MO diagrams for these oxides are given in Fig. ; 1 . 



Fig. 1(a). Valence molecular orbital energy level diagrams for Mg0 4 ~* and Si0 4 ~ 4 * 

(with j eigenvalues set equal). 

;Fig. 1 (b). Valence. molecular orbital energy level diagrams for MgOj l0 , A10 e -9 and 
SiO # ~ 8 t (with l< l3 eigenvalues set equal). 

In Tables l and 2 the eigenvalues and percentage atomic orbital composition of 
the molecular orbitals are given for the anions SiO.,~ and Si0 6 _1 *. The basis sets of 
atomic orbitals employed were the best atom minimum basis sets (BA MBS) of 
Cdementi and Raimondi (1963) and the M — 0 distances are given in Fig. 1. A 
tabular representation (Tables 1 and 2) is required in addition to the conventional 
MO diagrams given in Fig. 1 since the information on %AO composition is as irnpor- 
.tarit for our purposes as is the information on orbital eigenvalues. In the figures and 
•tables the MO's are labeled by symmetry type and by number within a symmetry 

* Distances used were R ( Mg O) -- 1-89 A. Ii ( Al - O) .« l-7-> A. W(Si- -O) = 1-60 A, 

f Distances used wore R(Mg— O) = 212 A. 7?(Al O) •- 1-93 A. /?(Si - -O) -- 1-80 A. 
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Table 1 . Eigenvalues ami atomic orbital percentage compositions of filled 


molecular orbitals for Si0 4 * 

Orbital* 

(eV) 

AO% composition t 

5/ g 

•29-2 

I % Si3p, 90 % 02/) 


27-9 

100% 02/) 

le 

26-2 

100% 02 p 

5a, 

23-3 

17% Si3*. 1 % 02.1, 82%02;o 

4t 2 

22-9 

l9% Si3/>, 7l%02/> 

3t t 

49 

3%Si3/). 97 % 02* 


3 1 

3 % Si3*. 97 % 02* 

2f 2 

-87-8 

99 % Si2/i, 1 % Si3p 

3a, 

-1371 

100 % Si 2* 


-508-2 

100% Ols 

2a i 

-508-2 

100% 01* 

lttj 

-1849-9 

100% Sil* 

* Orbitals 

are labeled by syinmetr 

v ty |K‘ and in order of increasing 

energy within 

each symmetry type. 


f The percentage of AO i in MO p 

is defined as 


X C pk S ik l2{ 

■ft . 


where the sum on k is over all AO's. 


Table 2. Eigenvalues and atomic orbital percentage compositions of filled 


tnoleculnr orbitals for Si() 4 "* 

Orbital 

(oV) 

AO % composition 

5*1.. 

53-5 

2% 02*. 97 % 02/> 


53-4 

100% 02/) 

^10 

53-3 

100% 02/, 


51-7 

100% 02/) 


49-1 

100% 02/) 

*hu 

47- 1 

14% Si3/>. t % 02*. 85% 02 p 

5°io 

46-2 

13 % Si 3*. 4 % 02*, 82 % 02 p 

-V . 

31-9 

100% 02* 

Shu 

30-4 

1 % Si3/). 99% 02* 


27-4 

•> % Si3*. 95% 02* 

%hn 

— 04-0 

98%Si2/). 2%Si3 p 


— 11 2 0 

99 % Si2* 

. 2a u 

-482-0 

100%0i* 

le„ 

-482-0 

100% Ol* 

K iu 

-482-0 

100%01* 

la t« 

-1827-6 

100% Sil* 


type in order of increasing energy. The symmetry types are determined by the irre- 
ducible representations of the tetrahedral or octahedral point group. Because of the 
excess negative charge on the systems, many of the filled MB's have positive eigen- 
values. However, we shall be interested only in differences of eigenvalues, so this 
point will cause no difficulty. 
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The valence M.() 4 s for Sh >., -4 can be divided into three groups: (l) -by, and 3 1 % 
which are -predominantly U2s orbitals, (2) 41., and ~>u , which are ( 7 -bonding orbitals 
with both Si and O character; ('■',) la, and 'it., which are non-bonding and weakly 
7 , --bonding predominantly (>2p orbitals. For Sit) K g , a similar division may be made; , 
with -b/j,,, .‘1/ l„ uiid' 2 rg being mainly 02 s orbitals: *>((,„ and 4h„ ( 7 -bonding orbitals; 
and the five highest orbitals non-bonding or weakly --bonding orbitals. More accu- 
rate SCT calculations on tetrahedral oxyanions (Gkucs, at id., 1071) indicate that 
the MO orderings in Table 1 are probably correct, except for the order of ~>t 2 arid \t v i . ' 

The relative energy of these orbitals may be in error by several eY. The orbital order X 
in Si0 6 ~ 8 should also be correct for the most part. 

Assignment of X-ray Emission Spectra " 

A detailed discussion of the assignments of the SiK„ K fi , L 23 anti QK, spectra 
for a-fjUJi'tz is given elsewhere (Tossec.e, 1072) . Here the assignments of the. MK* 
and K. and OK, spectra for J / -- Mg, A1 and Si in tetrahedral and octahedral co- 
ordination with oxygen are presented and calculated energies and intensities are 
compared with experimental results, where available. Group theory yields the 
selection rules needed to identify the allow transitions. 

MK a Spectra . ; 7 ; -g V 

The MK, transition, nominally M2j> — Mis in character, is identified with the 
ll 2 MO in tetrahedral and the 2t u MO in octahedral coordination. The K, line is 
found at higher energies in oxides than in the corresponding metal. The calculations ’ : . b 

give A'jjk.. values 2 to .'i eV larger in the oxides than in the free atoms. This trend ;|s . : 

in' the right direction but is substantially larger than the experimental metal to oxide 
. difference, which is about <>•."> oY. Experiment also indicates that E K increases with 
increasing M — 0 distance and with increasing coordination number of the metal, ; ' . : 

These trends are also reproduced by the calculations. Octahedral oxides are calcu- 
lated to have higher E MK . than tetrahedral oxides by about 1 eV, which is, however 
several times larger than the experimental difference, The higher value of E Kj in 
octahedral relative to tetrahedral oxides can be related to the higher Mulliken charge 
on the metal atom in octahedral coordination. In another place (Tosssbr., 1972) 
the K, sliifts for Si mid Al are discussed, in more detail and it is shown that the AIK, 
peak shift, as a function of /’(Al -O) for A1 in tetrahedral coordination (as measured . 
by and fb:iNO'.::v, M)72), is given by SCF XE.MU calculations to witliin j 

about a factor of 2 . | 

MK P Spectra ! 

In Fig. 2 spoct’-a are presented for and MgO. In the SiO., spectrum we 
observe three features; n satellite peak K r , a rna : n peak IC, and a shoulder on the 
high energy side of K,. In Mg' > the main Jv, peak is split into two components. :■ 

Alignments and calculated energy's and intensities for MK /( spectra are given in j 

Table 2. Intensities w ealeiilatod from the % 8i3j» character in the iudi vidu-d 310’i.'' J 

Maxxe (l'i)7i.') c ! !i:»!oy<Ml a similar nivtliod for liydrocarhou CK, spectra. Agreement i 

between calculation a id experiment is fairly good. For example, the relative energies 
and intensities of the three peaks in the K. /; spectrum of x-ipnvrtz are predicted with 




1640 1033 1830 1822 1830 ISIS 1810 I3C5 1300 1295 1290 1285 !280 

E. eV E * eV 

(a) - (b) 

Fig. 2. Experimental* MK^ X -ray emission s pectin for (a) x-rpmrt/ find (b) 

pcriclase. 


Table 3. Experiment;.! and er leuleted relative energies (i;V) and intensities in MK.i spectra 



■ Exp. • 


Otic. 


Ex,.. 


Calc. 










S*<\ 




AE 

/ 

As.-rig:ini«*iit A# 

/ 

AE 

/ 

Alignment 

AE 

/ 

K/r 

— 13'7, 
- 14*4 

13 

3(, 1«, -17-9 

III 

— 



-160 

9 


: 0 

100 

4l a -► l/i, 0 

too 

— 

; — . 

la l S 

0 

100 

wUmtldi.-r 

4 o*o 

13 

5/,— l.i, ~0*3 

5 

' — 

: — • 

*1. 

6-4 

1 




\l t 4 




AID, 




AE 

l 

A E 

i 

A £ 

i 


AE 

I 

K/»* 



-180 

17 

-10*8 

29 


— 17*1 

8 

Kp 



0 

100 

0 

100 


0 

100? 

shoulder 



4*7 

3 

— 

— 


4-7 

2 


AE 


Mg T d 




Mg O, 




I 

A E 

/ 

A E 

/ 


AE 

I 

Kir 



-18-0 

24 

— 10*2 

25 


-18*2 

28 

Kfi 



0 

100 

0 

100 


0 

100 

shoulder 



2-9 

2 

— 

— 


2-8 

0 


* Experimental values are taken from the tables of C«<:h (1970) and the figures of Dodo and Guess (1968). 


reasonable accuracy. The K r peak is assigned to the 3U MO of mainly 02? character, 
yielding a peak of the proper energy and intensity. The high energy shoulder on the 
main K t peak is assigned to the Ht., M(), yielding an acceptable energy and intensity. 
It is not necessary to assign this peak to a symmetry forbidden transition from an 
orbital of e symmetry, as in Dodd and Glexx (1908). The separation of main \\ B 
peak and the high energy shoulder is calculated to be 1 6eV smaller for tetrahedral 

'* From Dodd and Glenn, 10(58, 
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-Si than for tetrahedral A1 in cjp.se agreement with the experimental value of 1-5 eV 
(Brown ct ah, 1969). 

Ubqh (1970a, b) has used c/ualitative MO theory to explain the observed vuria- 
.tion in the K.-K.,. energy difFercnce (here designated A) within n series of oxides as 
..the cation changes from Mg to 1.-1. He observed that A decreased by 3-.7 eV along 
this series and that the M3p orbital also decreased in energy; i.e. became more 
strongly bound, with Chip more stable than Mg 3p by about 8 eV. The observed 
decrease in A was interpreted as a lowering in energy of the 4/ 2 MO of mostly M3p 
and M3s character relative to the 3/ 2 MO of mainly 02 s character, which decreased 
in energy at a slower rate. Since the variation in M:5p energy is about two or three 
tjmes as large as the observed A variation, the h 2 MO will jwobably be less than half 
;M3p in character. 

One also observes a change in A as the ligand, L, coordinating the central metal 
varies, with A increasing along the first row from 0 to F. The interpretation of Urch 
( 1970a, b) indicates that this occurs because of the variation of the energy of the L2$ 
and L2p orbitals. The difference of L2p and Lis energies for first row atoms increases 
at a rate of about 5 cV per unit increase in nuclear charge, similar in magnitude to 
the observed shift in A. Also, A values for oxides and fluorides range from 13-7-16-9 
and 19-21, respectively, in good agreement with the L2p-L2s energy differences of 
14-8 and 20 4 for 0 and F, respectively. This indicates that for oxides and lluoricles 
the 4 1, 2 MO (and its octahedral analog 4q„) is predominantly L2p in character, in 
agreement with the results presented in Tables l and 2. The 4/. 2 -3G energy difference 
is then predominantly a L2p-L2v difference, modified by M3p participation in the 
4< 2 MO. 

Since differences of L2p and L in energies obtained from a minimum basis set 
atomic SCF calculation are substantially larger than the experimental values tabu- 
lated by Urch ( 1970), good absolute values of A can not be expected from the calcula- 
tions. However, it is found that relative A values are given well. The quantitative 
MO results reported here give a decrease in A of eV r per unit increase in nuc- 

lear charge, in good agreement with the experimental value of about 0-6 -1 -0 eV/unit 
charge. However, the experimental data arc not good enough to establish whether 
the trend is as uniform as it is here calculated to be. 

The value of E MK<|/ is also found to decrease slightly from octahedral to tetra- 
hedral oxides in our calculations, in agreement with the results of Dodd and Glenn 
(1908). The calculated decrease is around 2 eV T for Mg, increasing to 3 eV for Si. 
This is several times larger than the experimental value. Also, calculated E ni - 
values are lower for tetrahedral oxides than for octahedral by about half the magni- 
tude of the K 0 ' shifts. The calculated SiK^ peak shift between SiO," 4 and SiO ( .'• ,, is 
.calculated to be about 40 per cent larger than the experimental difference between 
(4-quartz and stishovite (White and Gibbs, 1067). Our calculations indicate that A 
is lower in octahedral than in tetrahedral oxides by about l c V. Experimental data 
are presently not available to test this result except for those oxides in which the 
splitting of the main K t peak renders such analysis difficult. 

A rather curious feature of the MK^ spectrum is the splitting of the main peak 
by several electron volts which occurs in several oxides. A related observation is 
that peak widths are larger in octahedral oxides in which such splitting is more 



J. A. To.ssku. 


commonly observed. Tbe peak splitting in eonindnm (x-ALOj), with Al in octa- 
hedral coordination, has been interpreted by Dodd and Clucxx ( 1008) as a result of 
transitions from the 4f, tl MO of Al 3pa character and t he 51 lu MO of Al ?>p- character. 
Although our calculated 5t lll -4l lll energy differences could be reconciled wit!) this 
interpretation wo must reject it since the %M character in the ~>l in MO is very 
small and the intensity of the peak generated by this MO is therefore quite low. The 
5t lu MO in octahedral symmetry is chemically almost identical to the 51 s MO in 
tetrahedral symmetry which has been assigned to the weak high energy shoulder in 
the a-quartz spectrum. Another possible explanation for the splitting is a dis- 
tortion of the A10 fi octahedron. However, SOF XEMO calculations show that even 
a distortion of 0-- A would yield a Kj splitting of less than l eV. The explanation 
must then lie in the interaction of neighbouring units in the structures. There is an 
important difference between the type- of linkage of the fundamental metal-oxygen 
units in a-quartz and, for example, pcriclase (Mg 3) which shows a large splitting. 
In a-quartz the Si0 4 units are linked by corner sharing. In periclase the linkage is 
by edge sharing. Interaction between Si atoms in x-quartz is therefore small since an 
oxygen core occupies the region between t he two silicons, In MgO the oxygen ligands 
are at some distance from the Mg -Mg intemuclear direct ion. The interaction of the 
Mg atoms is therefore essentially unimpeded. This geometrical difference is shown 
in Fig. 3 for the dimer geometries used in the calculations to described below. 



Fig. 3. Geometries employed for dimeric units Su0 7 -8 nnd Mg 2 O 10 18 . 


To test the qualitative hypothesis presented above SOF NEMO calculations were 
performed on the dimeric units Si./) ; 0 and Mg./) The distances and angles 
chosen yielded essentially identical M — M distances [7?(8i -Si) — 3-02 A, i?(Mg— - -Mg) 
= 3-00 A]. The Mg valence orbitals were only slightly more diffuse than those of Si 
and yielded overlap integrals larger by only about 10 per cent. On the other hand, 
an extremely large difference was observed in M — M overlap populations [designated 
n( A — B)], with ri(8i — Si) •— — 0*U6 and »(Mg -Mg) —0-4 f, indicating qualita- 
tively different M M interactions in the a-quartz and periclase models. In the 

periclase model the M — M interaction is very strongly antibouding. 

The ‘44a like' MO’s in Si/ > 7 " are also spread over a smaller energy range than are 
the ‘4/ lu like’ MO’s in Mg/.)^ 10 . However this may be clue partly to the smaller 
number of oxygens engaged in bridge bonding, in the Si clitncr. The most important 
result of the calculations is clearly the difference in M— -M overlap populations. 
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Although our dimeric results can only give a crude picture of the electronic 
structure of a complex solid, they quantitatively support the hypothesis proposed. 
The M — M interaction in edge shared Al'g.f is found to be inuch greater than that 
in corner shared Si 2 v) 7 The large metal-metal interaction .in MgO is hypothe- 
sized to be the cause of the MgK^ splitting. 

OK a 

Our assignment of the OK, emission spectrum of x-quartz (O’Nions and Smith, 
- i 1)7-1 ; Smith and O’Nions, 1972a) is given in Table 4. The 531 cV peak in their spec- 
trum has recently been interpreted by them as mainly a reflectivity spike, although 


Tabic 4. Experimental and calculated relative energies (eV) and intensities 
of OK a spectra for a-qunrtz 


Exper. 


Calc. 



A E 

/ 

Valence orbital 
assignment 

AA 

I 

0 \ 


5t, 

0 

99 

-10 

4 

1*1 

-1-2 

100 

-2-8) 


le 

-30 

67 

— 6-8} 

I 

5a, 

-5-9 

27 



4 *2 

-6-3 

71 


some characteristic intensity remains after correction for this effect. The remaining 
characteristic intensity probably arises from a double ionization process, which is 
outside the scope of our calculations (Fischer, 1965). If this peak is disregarded, 
agreement between calculation and experiment is quite good. The division of the 
spectrum into a high energy and a low energy band is reproduced with reasonable 
energies and intensities. 

The high energy OK, band, at about 525 eV, arises from the non-bonding MO’s 
and lies slightly higher in energy in Si oxides than in A1 or Mg oxides. Also the lower 
energy peak, (520-522 eV), arising from the a bonding orbitals, is separated from the 
high energy peak by different amounts in the Si, A1 and Mg oxides (Smith and O’Nions, 
1972). The separation increases along the series from Mg to Si, consistent with the 
calculated MO eigenvalues presented in Fig. 1. Experimental values of this separa- 
tion (Smith and O’Nions, 1972) are 2-9 eV for pcriclase, 4-1 eV for corundum and 
4-9 eV for a-quartz. The calculated separations between the (lfj, lf ly ) MO and the 
(4t 2 , 4 t lu ) MO for MggOi, 10 , A10 6 ~ 9 and Si0 4 " 4 are about 2-8, 4-6 and 5-1 eV, respec- 
tively. The agreement of calculation and experiment is therefore rather good. The 
increase in the peak separation along the series Mg to Si is indicative of an increase 
in the stability^ of the a bonding orbitals along the same series. The analysis above 
is clearly dependent upon the decision to ignore the anomalous 531 eV peak and its 
success supports this decision. 

Electronic Structure ok Mg,. A! and Si Ox van ions 

Having established the ability of Sf.'F XRMO calculations to reproduce trends 
in X-ray emissions spectra we can have some confidence in the other properties 
.obtained from them. Here some features of the calculated electronic structures are 
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discussed briefly. The variat ions of several charge distributions parameters for the 
systems considered are presented in Table 5. The metal charges are found to increase 
uniformly along the series Mg, Al, Si and in passing from tetrahedral to octahedral 
coordination. The variation in the M3s atomic orbital population is fairly small and 
the total number of electrons on the metal is fairly constant. However, a uniform 
increase in is observed, indicating that the ratio of K f intensity to K a intensity 
should increase along the scries Mg to Si. The charge on oxygen also decreases 
uniformly along the same series. However, this may well be a result of the decreasing 
charge on the oxyanion. The increase in unbalanced electron repulsions as the an- 
ionic charge increases will also tend to place too many electrons on the ‘outside’ of 
the molecule, i.e. upon the oxygens. If this effect is real it would manifest itself in 
higher energies for t he OK, peaks in oxides of Si than in the corresponding oxides of 
Al or Mg. As mentioned before, the experimental data is in general agreement with 
this result (Smith and O’NTons, 1972). However, the relative oxygen charges must 
be uncertain until larger clusters are studied. 



Table 5. Charge 

distribution 

in MO| 

and MOg anions 




Mg0 4 -« 

MgOj 10 

aio 4 -' 

- AKV 8 

SKV* 

sio 6 - 8 


106 

1-26 

1-88 

1 -90 

2-58 

2-83 

<?(0) 

— 1-77 

-1-88 

-1-72 

— 1-83 

-1-64 

-1-80 

^M3s 

0-21 

0-22 

0-20 

0-25 

0-39 

0-37 


0-83 

0-67 • 

100 

0-93 

IOC 

0-87 


14 5 

112 

1G-7 

133 

18-8 

13-8 


* The numbers given tiro %M3j» in cither the or the 4<j u MO, depending upon the coordinn. 
tion symmetry. ' 


The higher values of E MK(( in octahedral oxides, as compared with tetrahedral, 
seems to be connected with their lower covalency. The amount of M3p per cent 
character in the 4f lt< MO in octahedral oxides is consistently smaller than the M3p 
per cent character in the analogous tetrahedral point group orbital,. the it v The 
calculations also show a gradual increase in metal orbital participation in the a 
bonding orbitals along the series Mg to Si for both tetrahedral and octahedral oxides 
(see Table ~>). Figure l shows a substantial spreading of the MO energy level diagram 
along the same series. Noniterative NEMO calculations using neutral atom a’s show 
the same trend. Therefore both MO compositions and MO energies indicate a sub- 
stantial increase in covalencv along the series Mg to Si for both tetrahedral and 
octahedral oxides. Figure la shows that as the covalency, increases both the it 2 and 
3 1 2 MO’s (as well as their octahedral analogs shown in Fig. lb) are lowered in energy 
relative to the l/j. However, the it.,, the more strongly a-bonding of the two, is 
lowered by a larger amount and the value of A, the K f -YL f . separation, is therefore 
decreased. The value of A is then a fairly direct measure of covalency. 

Conclusion 

The quantitative MO calculations presented here have been shown to be capable 
of assigning the M/v ll 2 peak and all the peaks normally observed in the MK # and OK a 
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-spectra without invoking symmetry forbidden transitions or making ad lac assump- 
tions about bonding type. Trends in energy differences fesiilting from varia- 

tion of metal, ligand and geometry have been reproduced by the c alculations and 
Urcii's (1070a, lDTob) interpretation of MK ;; spectra, which is qualitatively correct, 
has been made more quantitative. The splitting of the main MK^ peak in several 
oxides has been related to the type of- linkage of t-heir fundamental metal-oxygen 
units. This demonstration of good agreement between calculation and experiment 
for X-ray emission spectra suggests that the calculations reported here will generally 
yield an accurate picture of t he electronic structure of oxalic compounds of Mg, .A 1 
and Si. The calculations show the gross MO structure "of the various oxides to be 
quite similar. In all cases the valence Ml Vs can he divided into three groups, defined 
in terms of energy and composition. Quantitative differences in the calculated 
electronic structures correspond closely with X-ray emission spectral results. In 
particular, A is shown to be closely related to the covalency, which increases; sub- 
stantially along the series Mg to Si. 
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molecular orbital (mo) interpretation OF THE METAL Ko, metal 
L AND OK X-RAY EMISSION SPECTRA OF SOME OXIDE AND SILICATE 
MINERALS 

Tossell, John A. and Vaughan, David J., Department of 
Earth and Planetary Sciences, Massachusetts Institute 
of Technology, Cambridge, Mass. 02139 
Molecular orbital calculations, using an approximate SCF 
LCAOMO method, are presented for several oxides of Mg, Al, 

Si and Fe. The valence MO ' s for the tetrahedral third row 
oxides can be divided into three groups: (1) 4a^ and 3,t 2 ; 

02 s orbitals, ( 2 ) 4t 2 and Sa^ o-bonding orbitals, (3) le, 
ltj and 5 1 2 ; non-bonding 02p orbitals. A similar division 
can be made for the octahedral oxides. The quantitative MO 
diagrams are used to assign the observed MKg and 0K a X-ray 
emission spectral peaks for several oxide minerals.. The 
trend in the MKg-MKgi peak separation is calculated for the 
series MgO, AI 2 O 3 , Si0 2 arid results are in agreement with 
experiment. The peak separations in the 0K Q spectra show an 
increase along the same series, also given accurately by the 
calculations (exp. 2. 9, 4.1 and 4.9ev, calc. 2. 8 , 4 . 6 and 
5 . lev). The MKg and 0K a peak separations are shown to be 
quantitatively related to the degree of metal-oxygen coval- 
ency. Fe L X-ray emission spectra are presented for several 
silicate minerals containing iron in different oxidation ' 
states and coordinations. The main L Q band is composed of 
four peaks which are assigned to the four valence MO's of e 
and t; type which possess substantial amounts of Fe3d char- 
acter. An additional peak at lower energy is assigned to a 
MO which is predominantly 02s but has 1-21 3d character.’ 
This peak is. found to be significantly higher in energy in 
ferric compared to ferrous compounds. 




X-ray Photoelectron, X-ray Emission and UV Spectra of 
SiO ? Calculated by trie SCF X o< Scattered V/ave Method 


J. A. fossell and D. J, Vaughan 
Dept, of Earth and Planetary Sciences, 

Mass. Inst, of Tech. f^/f 

and . 

<r K. H. Johnson 

Dept, of Metallurgy and Materials Scier.ce, 
Mass. Inst - , of Tech. /— -h^.11 

(for submission to Chem. Phyo. Letters) 



Abstract 

Photoionization ; X-ray emission and UV excitation 
energies are calculated for Si0 2 by the SCF Xc* method and the 
transition-state procedure. In all cases agreement between 
calculation and experiment is good. The SiO^”^ cluster 
is found to be adequate for describing localized excitations 
in quartz. The valence orbitals of this cluster are found to 
be naturally separable into three sets (1) 02s, nonbonding; 
(2) Si3s,3p and C2p, sigma bonding; (3) 02p nonbonding. The 
participation of Si3d type functions in the high energy 02p 
nonbonding orbitals is found to be very small . 



Research by Aviva Brecher 

An experimental study of magnetic properties of lunar rock 
(breccias) chips has been undertaken, complementary to current 
studies on them by Prof. Burns and coworkers. To date, I have 
examined the Apollo 16 samples 60315.51 and 62295.27. My non- 
destructive magnetic measurements are aimed at: 1) determining 

the intensity and stability of natural remanence in lunar rocks 
under study; 2) assessing the weight fraction of magnetic phases 
and the level of magnetic saturation; 3) identifying the minerals 
(oxides or native iron) carrying the remanence; 4) determining 
the size spectrum and associated magnetic behavior of carrier 
grains as a function of temperature. These have been achieved 
partly by studying coercivity, saturation moments and reduced 
saturation ratios. 

Currently this experimental work entails great unnecessary 
effort, being carried out on improvised equipment, heavily used 
by workers in other departments at M.I.T. It is hoped to develop 
soon an adequate laboratory facility for studies of lunar and 
meteoritic paleomagnetism . 

Theoretically, I intend to reexamine evidence on magnetism 
of lunar crystalline rocks, which is not consistent with the cur- 
rently favored interpretation as thermoremanence , acquired in a 
global lunar field of ^ 2000 y. This interpretation has far-reaching 
implications regarding the thermal history of the moon and based 
only on remanence of breccias. Experimentally, it is hoped to 
carry out ferromagnetic resonance studies of lunar rock in order 



to determine "in situ" relative volume fractions of grains of a 
given elongation. This will test Fuller's suggestions that; 
remanence in lunar rocks and soils, although blocked at low 
temperatures in fine iron grains is stable under AF-cleaning be- 
cause of their shape anisotropy. 

It is also hoped to study the acquisition of inverse thermo 
remanence in lunar rocks, and the effect of thermal cycling 
(-80 to +140°C) , stress and shock on it. 
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D. Leavy 

Progress report on: Lunar Electrical Conductivity Project 

The aim of this study is to invert the magnetic field spectrum 
( £ < 5 x 10 2 hertz) observed in the solar wind and on the lunar 
surface in order to deduce the internal electrical conductivity 
of the moon. The electrical environment of the moon is modelled 
by a perfect conductor on its sunlit side and an empty cylinder 
on the downstream side. (See Fig. I.) A numerical scheme, de- 
tailed in App. I, is adopted in order to calculate the field induced 
by a layered moon. 

Fig. 2. compared the data obtained by Sonett, et al . with 
the Apollo XII magnetometer with the field calculated by Sill and 
us for the same model. We see that as the frequency increases. 
Sill's calculation becomes inaccurate so that a large part of the 
spectrum cannot be used to model the moon response adequately. 

In an attempt to model the asymmetry in the lunar plasma 
environment, Sonett et_ al. have recently modelled the moon as an 
infinite conducting core embedded in a non-conducting crust. We 
have tried to reproduce their results [Fig. 3, 4] within our theory 
and the result seems to agree quite well. Of course, our analysis 
is not constrained to such a model so we might expect that it will 
be a better tool to find the lunar conductivity. 

Presently, we are trying to develop an inversion scheme and a 
resolution analysis of the data. This work should be completed, 
we hope, in early 1973. 
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Appendix I 
Mathematical Theory 

We consider a magnetic fluctuation of wavenumbers 



impinging on the moon normally to the axis of the cylinder. 
X « e • , 
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The field inside the cylinder is given as the sum of two parts: 

1) An inhomogeneous field satisfying the continuity of 
the tangential electric field at the plasma-cylinder 
boundary; 

2) A field having no tangential electric field at the 

solar wind-cylinder boundary, but satisfying the continuity 
of the tangential electric and magnetic field at the 
moon-cylinder boundary. / 

Formally we have 
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where and X 2 ^ are respectively the root of the cylindrical 
Bessel function and its derivative 

X - iki 

an'd o and 3 are the coefficients of the TM and TE modes respective- 
n n 

iy. 

Inside the moon, the field is expanded in a sum of toroidal 
and poloidal mode. 

At the surface of the moon, we must meet the continuity of 
the tangential electric field on the front side. The continuity of 
the tangential electric and magnetic, field must also be met inside 
the moon and on the back side. 



Imposing these conditions at a certain number of points 
around the surface of the moon and including in the resulting 
system of equations only a finite number of modes, we can 
solve for the field by a least-square technique. The number Of 
modes included in the analysis is determined by the criterion 
that only a few percent mismatch in the boundary condition is 
allowed. 
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PHYSICAL SCIENCES 

On the Nature of the Tycho 
X-ray Source 

Of the five X-ray sources which are identified with remnants 
•of supernovae, the Tycho 1572 source is of particular interest. 
In this case the X-ray emission cannot be considered as the 
continuation of the synchrotron radio spectrum because it can 
easily be shoes n that, independent of the mechanism of injection, 
the spectrum must be "bent” in t he frequency region s'< 10 13 
Hz. Also, it cannot be assumed that the observed X-radiation 
is thermal and caused by hot plasma behind the front of shock- 
wave produced by the expanding envelope of a remnant. 
This mechanism qualitatively explains X-radiation from 
Cassiopeia A if the undisturbed density n 0 of interstellar 
matter in the vicinity of the explosion is about 30 cur 1 . 
This value of n 0 is in excellent agreement with the density of 
the HII region which is observed around Cas A (ref. |). But 
in the case of the Tycho remnant one must assume that 
n 0 ~ 5 cm -3 if the radiation is thermal. Such a value of n 0 
contradicts the observed large dimensions of this remnant, 
which, according to all indications, has exploded in a relatively 
rare medium between the clouds of interstellar matter, where 
n o <0.\ cm ~ 3 . This excludes the possibility of explaining 
X-radiation from Tycho as thermal radiation from hot plasma. 

It seems probable that this source should be identified with 
a young pulsar which originated after the explosion of the 
Tycho supernova. Because of "unfavourable” orientation of 
the rotation axis of this pulsar, coherent radio emission is not 
observed. 

But there are good reasons to believe that young pulsars 
should generate synchrotron radiation in the optic and X-ray 
ranges (similar to NP 0532, ref. 2). The power of this radiation 
must be a very steeply diminishing function of the pulsar’s 
age. It is essential to note that the intensity of synchrotron 
radiation of NP 0532 between pulses is not zero. For example, 
according to Vishvanathan 3 , the intensity of radiation between 
the main and secondary "peaks” is about 1.5% of the maximum 
intensity in the main ‘’peak”.. According to Bradt et a!.*, in 
the X-ray range the “intcrpulse” intensity is near 10%. In 
the 100-400 keV range this intensity is sufficiently larger 
(I. Kurfcss, unpublished work). Consequently, one should 
remember that there is a “quasi-isotropic” component in the 
X-radiation of young pulsars. The time-averaged power of 
this component is at least equal to the “pulse” component. 
The X -ray flux from the Tycho supernova remnant is -25 
times smaller than that from the Crab Nebula 5 and, conse- 
quently, it is half the X-ray flux from NP 0532. Since the 
distance to SN Tycho is approximately 2.5 times that to the 
Crab, the observed X-radiation from the Tycho remnant can 
be explained if the power from the expected pulsar, which is 
situated inside Tycho, is three times that from NP 0532. This 
can be explained naturally by the youth of Tycho’s pulsar. 
Recently Henrikscn et at. (unpublished. work) suggested that 
X-radiation from Tycho and some other sources can be inter- 
preted as thermal radiation of his crust heated due to dissipa- 
tion of rotation energy of a neutron star, in our case it should 
be supposed that the crust temperature is 25 x10*’ K: But in 


this case the linear dimensions of the source must be unaccept- 
ably small— above I km. 

The observable tests of our hypothesis arc as follows: 
(a) measurement of angular dimensions of -Tycho’s X-ray- 
source — according to our hypothesis, the upper limit of angular 
dimensions of this source must be essentially smaller than 
those of the radio source; (h) discovery of a small periodic 
modulation of X-ray flux from Tycho caused by expected 
anisotropy of radiation of Tycho’s pulsar; (c) a study of l.Slh 
to 19th magnitude stars in the central part of the Tycho radio 
source wit!; the object of discovering an optical pulsar with 
periodically modulated flux; U/) discovery of a "point”. X-ray 
source in the remnant of Kepler’s Supernovae 1604. 

The “quasi-isotropic” component of pulsars’ synchrotron 
radiation must play an essential role in the physics of the Grab 
Nebula. Evidently, this radiation generated by relativistic 
electrons closed in the radiation belts of NP 0532. The 
spectrum of this radiation is essentially harder than in the 
main pulse. This can explain the expected relation between 
puisars and nebular synchrotronous spectra 2 . 

It should be noted that in connexion with Pacini's explana- 
tion of “glitch" in NP 0532, the systematic observation of 
"quasi-isotropic” component synchrotronous radiation of 
NP 0532 is desirable. 

I. S. Shklovsky 

Institute for Cosmic Research , 

Moscow 

Received February 23, 1972. 
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Lunar Electrical Conductivity 

Sonett et a/. 1,2 have inverted sunward lunar surface mag- 
netometer data to obtain a lunar electrical conductivity profile, 
using the poloidai magnetic response to a spatially uniform 
source field of a radially structured conducting sphere immersed 
in an infinitely conducting plasma 3 - -1 . The experimentally 
determined lunar response expressed as the tangential mag- 
netic spectral amplification factor is shown in Fig. 1 together 
with the best fit poloidai responses for a two-layer and a 
multi-layer. Moon and their corresponding electrical con- 
ductivity profiles 5 . • 

The spatial symmetry imposed on the boundary conditions 
by this model of the interplanetary magnetic field-Moon 
interaction neglects the solar wind plasma-Moon-void 
geometry and the finite applica .ion time of the interplanetary 
magnetic field which result from the solar wind flow past the 
Moon. We have studied qualitatively the effects of including 
these boundary condition asymmetries through numerical 
examinations of an asymmetric cylindrical model of the inter- 
action (Fig. 2). In our model a cold, infinitely conducting 
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Vil( , ( j plasma is incident with velocity v 0 x and the solar 
' j l,,,nc magnetic field II is aligned along the axis z ot a 
- JocjI Moon h . 

' j shows the numerically determined sunward surface . 

1 w of the asymmetric model of the interaction for a 
—rvru't conductivity model with a highly resistive thin 
" ff'b'cr added to eliminate the toroidal mode. The anaiytt- 
■ v ', determined cylindrical analogue of the poloidal response 
s. runic trie model 0 is also shown. Comparison reveals.. 



fc. t Experimental tangential amplification factor from 
v -.cit ,7 «/.\ shown with the calculated poloidal responses lor 
re inset conductivity profiles reproduced in part from ref. 5. 



***- ’ Asymmetric model of solar wind borne magnetic field 
interaction with cylindrical Moon. 



first, the depressed asymmetric response at low frequencies 
resulting from diffusion of field lines around the core into the 
anti-solar cavity on a time scale poOr.Lv where C~ R M is 
the scale length of the diffusion process; and second, the 
flattened asymmetric response for high frequencies .which 
results from application of the incident field to the sunward 
lunar surface on a time scale of R M /r o- 
These features of the response of the asymmetric model are 
also present in the experimentally determined response (Fig. 1). 
An attempt to fit these features using the poloidal response of 
the symmetric model of the interaction resulted in the compli- 
cated conductivity profile obtained by Sonctt et at. We also 
note that the response of the asymmetric model is insensitive 
to increases in the core conductivity cii, from which vve con- 
clude that the experimental sunward tangential magnetic 
response docs not exclude a more conductive lunar interior. 

A. Colbert Reisz 
— Donald L. Paul 

Theodore R. Madden 

Department of Earth and Planetary Sciences , 

Massachusetts Institute of Technology, 

Cambridge, Massachusetts 02139 

Received February 9,. 1972. 

1 Sonett, C P., Colburn, D. S., Dyal, P.. Smith, B. F., Schubert. G., 

and Schwartz. K.. Sat are, 230, 359 (1971). 

2 Sonett, C. P., Schubert, G.. Smith, B. F., Schwartz. Kv. and 

Colburn, D. S., 1‘roc. Apollo 12 Lu/t. Sci. Conf (MIT, I '->7 1 ). 

3 Blank, J. L„ and Sill, W. R„J. Ceoplns. Kes., 74, 736 (1969). 

4 Sill. W. R., and Bi.tnk. J. L., J. Geoplivs. Res., 75, 201 (1970). 

5 Sill, \V. R., Bclkonun TM-71-2015-9 (1971). 

6 Reisz, A. C., Paul, D. L., and Madden, T. R., The Moon (in the 

press). • ,H' 


Lunar Electrical Conductivity— Reply 

Reisz et at. 1 consider, in the preceding letter, tvyo concep- 
tually important points regarding the lunar transfer function 
inversions which vve reported 2 ’*. That work used a sym- 
metric confining current theory so that the dark hemisphere 
diamagnetic cavity elVect was ignored, as vve pointed out 
at that time; Blank and Sill 4 estimated that the effect upon 
the sunward hemisphere transfer function is small. Recently 
Schubert and Lee (private communication) have computed 
the fields associated with a magnetostatic configuration where 
confining currents exist sunward of the Moon, upon the 
surface of a perfectly conducting core and upon the walls 
of the diamagnetic cavity. Their results are applicable here 
by relating their core radius to an effective skin depth at a 
specified frequency. Approximate damping depths are given 
in ref. 3. For frequency /<0.001 Hz and R c ~0 . 6 R m where 
/?„ is the core radius and R m the lunar radius, the cavity 
correction to front side data not adjacent to the terminator is 
about 4%. For / >0.01 Hz the correction increases to about 
20%. The effect upon conductivity is trivial for Re ~ 0.6 R m 
and is a factor of about 3 for/>0.0l Hz; the resultant change 
in estimated deep temperature is insignificant. For the shell, 
dominated by the higher frequencies, a combined change in 
shell depth and conductivity takes place but is beyond the 
scope of this letter. So we do not agree that field line diffusion 
through the crust and past the core can account for the very 
large difference in transfer function in Fig. ,3 of ref. I for 
0.002 </< 0.006. VVc note also that the skin depth at /re 0.002 
Hz for the core conductivity* indicates a core radius, of only 
900 km, hardly sufficient to produce the effects described. 
Lastly, the difference between the two cases shown in Fig. 3 
of ref. 1 is non-monotonic, vvhich does not seem likely i'or 
various reasons. For example, higher harmonic effects ((> 1) 
are unimportant below /~ 0.008 Hz (Fig. I). On the basis of 
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MORPHOLOGY OF SELECTED LUNAR FEATURES 
USING EARTH -BAS ED RADAR 

G.H. Pettengill and Pierre Biraben 

One of the major difficulties lying in the way of 
using the full resolution possible in the preparation of 
lunar maps by radar at the Harstack Observatory has been 
inaccuracy in the working ephemeris. While the intrinsic 
lateral surface resolution is on the order of 2 km with an 
associated height measurement accuracy of about 200 m, 
inaccuracies in the standard values for the lunar orbital 
motion, physical libration and rotation led to lateral 
discontinuities of up to 15 km when maps prepared at 
different times were compared. The corresponding dis- 
crepancies in surface height were of the order of 2 km. 

As a result, the assembly of composite maps of topography 
and radar scattering characteristics were impaired. 

In the course of other research at MIT, I. I. Shapiro 
and his students have developed a versatile planetary 
ephemeris program (PEP) which has some capability to accept 
observations of the Moon and to yield improved estimates of 
the necessary values listed above. During the course of the 
present reporting period, Pierre Biraben has developed a 
way of using PEP to produce corrections to the original 
working ephemeris which yield predictions of lateral 



positions to an accuracy of about 2 km. He has developed 
the necessary intermediate digital programming to use the 
outputs obtained from PEP to calculate selenocentric lati- 
tude and longitude of the lunar subradar point, the position 
and doppler angles, and the apparent rotational velocity of 
the Moon at 5-minute intervals throughout each observation. 

A total of 50 days of observation were included, extending 
from October 1970 to August 1972. Comparison of these 
improved predictions with the measurements shows a lateral 
consistency of about 1 km, with only one case approaching 
a discrepancy of 5 km. Topography of adjacent areas seems 
to be consistent to within 100 m. 

Further work under this grant has ceased as of 31 August 
1972 because of the transfer of the personnel to other 
research. 



Lunar Research Statement 
H. F. Hinteregger 

I have spent most of my time since the beginning of 
November 1972 constructing and testing 5 pairs of Dif- 
ferential Doppler Receivers (DDRs). These are semi- 
production versions of the 2 original prototypes which 
were developed by C. Counselman and myself in the preceding 
several months and which were successfully used in trial 
observations between STDN (NASA) Tracking Stations in 
Florida and California on Oct. 27 and 28, 1972. Three pairs 
of the new DDRs will be used within the next three weeks at 
the tracking stations in Madrid, Florida, and Ascension 
Island. In March, intensive observations using at least 
6 tracking stations will be started and I will be respon- 
sible for readying at least a sufficient number of DDRs 
for this purpose. 
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Lunar Research Statement 
R. W. King 


My research work since Nov. 1 has been related to the 
determination of the lunar physical libration from differential 
VLBI tracking of the ALSEP transmitters . Most of my time during 
this period was spent reducing the data taken in a two-day test 
performed Oct. 27 and 28 using the STDN stations at Goldstone, 
Calif., and Merritt Island, Pla. From these data we have ob- 
tained post-fit residuals at the level of about 50 cm in equivalent 
distance on the moon. These results, together with some analysis 
of anticipated accuracies in the determination of lunar baselines 
and the physical libration, have been prepared for presentation 
at the Conference on Lunar Dynamics to be held at the Lunar Science 
Institute in Houston, Texas, Jan. 14-17, 1973. 

In addition to analyzing the test data, I have spent some 
time in programming tasks related to the eventual processing 
of a significant amount of these data. In particular, the 
equations for the numerical integration of the 

rotation of the moon about its center of mass and the partial 
derivatives with respect to libration parameters have been coded, 
though not completely checked out. Finally, I have obtained from 
the Space Science Data Center at GSFC all of the lunar laser 
ranging data which lave been released by the NASA LURE team, and 
I have compiled these data in a format which can be used with 
our program. 

BK/jbu 
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Jack Pines 

i ; • . . - 

An ultrasonic seismic modeling set-up has been 
constructed and is in thej process of being refined. Both 
surface and body waves have been observed before and after 
the introduction of various configurations of scatterers 
in the medium to determine the effect produced by them. 

The indication at this point is that a small number of 
surface scatterers will destroy surface wave propagation while 
a large number of randomly placed scatterers will generate 
energy envelopes and dispersion characteristics similar to 
those observed in lunar seismograms. 



Gilbert Wayne Ullrich 

During the past four months, work was accomplished 
Comparing meteor impacts on the lunar surface to the 
impacts of man-made objects. Specifically addressed, was 
the question of the relationship between the seismic coupling 
coefficients of the two types of events. This information 
is important for it would allow the man-made impacts to serve 
as calibration events to determine, through seismic measure- 
ments, the energy of meteors impacting the surface. It 
was learned that, because of the lower peak shock pressures 
of the man-made impacts, ;they would be expected to have a 
higher seismic coupling coefficient. A more detailed study 
of the mechanics of impact events is required to estimate 
the relationship between the coupling of the two types of 
events . . 
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ABSTRACT 

Possible models for the thermal evolution of the Moon 
are constrained by a wide assortment of lunar data. In 
this work, theoretical lunar temperature models are computed 
taking into account different initial conditions to repre- 
sent possible accretion models and various abundances of 
heat sources to correspond to different compositions. Dif- 
ferentiation and convection are simulated in the numerical 
computational scheme. ri; 

Models of the thermal evolution of the Moon that fit 
the chronology of igneous activity on the lunar surface, 
the stress history of the lunar lithosphere implied by the 
presence of mascons , and the surface concentrations of 
radioactive elements, involve extensive differentiation early 
.in -lunar history. This differentiation may be the result 
of rapid accretion and large-scale melting or of primary 
chemical layering during accretion. Differences in present- 
day temperatures for these two possibilities are significant 
only in the inner 1000 km of the Moon and are not resolvable 
with presently available data. 

If the Apollo 15 heat flow is a representative value, 
the average uranium concentration in the moon is 65 ± 15 ppb. 
This is consistent with achondritic bulk composition (be- 
tween howardites and .eucrites ) for the Moon. 
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1. INTRODUCTION 

The development of the thermal regime is the most 
important factor in the evolution of a planetary body. 

For the Moon, as well as other planetary objects, magma- 
tization, differentiation, and other processes of evolution 
are controlled by the thermal history. The interpretation 
of geological, geochemical and geophysical data from the 
Moon requires temperature information since bulk properties 
of materials are dependent on temperature. The purpose 
of this paper is to investigate the thermal evolution and 
present-day temperature models for the Moon as currently 
constrained by the vast amount of data from the lunar 
missions . 

The most serious studies of the thermal evolution and 
heat balance of the Moon started with the pioneering efforts 
of Urey (1951, 1952, 1955, 1957, 1962). Other studies that 
followed (MacDonald, 1959, 1963a;Levin, 1962; Anderson and 
Phinney, 1967; Fricker et al . , 1967; McConnell et al., 1967; 
Hanks and Anderson, 1969) benefitted from and built upon 
those of Urey. The successful missions to the Moon have 
provided a large amount of new data and strong constraints 
on the present-day and past temperatures in the lunar interior. 
These data have led to renewed interest in evolutionary and 
thermal models of the Moon (Urey and MacDonald, 1971; 
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Papanastassiou and Wasserburg, 1971a; Hays, 1972; Wood, 

1972; Reynolds et al_. , 1972; McConnell and Gast, 1972; 

Toksttz et al. , 1972c; Tozer, 1972; Hanks and Anderson, 

1972). A common theme for most of these models has been 
to explain the early magmatism indicated by the lunar 
chronology and feldspathic crust, and the relatively cool 
present-day internal temperatures as suggested by the 
low electrical conductivity, high viscosity and paucity 
of tectonic activity. 

In this paper, we present lunar thermal models based 
on all data (including Apollo 16) and realistic computational 
schemes. The constraints which are most important for 
evolution models are discussed in Section 2. The compu- 
tational procedure and input parameters are given in Section 
3. In Section 4, the thermal models are presented, and 
this is followed by a discussion of the results. 
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2. CONSTRAINTS ON THERMAL EVOLUTION MODELS 

The data from lunar missions, especially those from 
Apollo experiments and samples , impose a wide range of 
constraints. Some of these (such as the electrical conduc- 
tivity, heat flow, seismic velocities and tectonism) 
pertain to present-day temperatures inside the Moon. Some 
others (such as the history of differentiation and large 
scale magma t ism, the existence of lunar mascons) indirectly 
point to the thermal state early in lunar history. 

These data impose some boundary conditions, and at times, 

what might appear to be conflicting requirements on the past 

and present-day temperature models of the lunar interior. 

For example, interpretation of electrical conductivity appears 

to suggest relatively low temperatures (Sonnet et al. , 1971, 

1972; Dyal and Parkin, 1971; Dyal et al. , 1972), in agreement 

with low seismicity (Latham et al. , 1971, 1972a,b) and high 

viscosity (Baldwin* 1971; Arkani-Hamed, 1972) ... In contrast 

to such information is the high value of surface heat flow, 

2 

33 erg/cm -sec (Langseth et al. , 1972). To better understand 
the limits of these major constraints, we discuss each one 
briefly, starting from those that relate to the initial 
conditions and going on to those that bear on the present- 
day temperatures. 


2.1 Lunar Differentiation and Choonology 

The Moon went through a period of major igneous activity 
and differentiation at an early age, during the first 1.5 
billion years after formation. A feldspar-rich, 60-km thick 
crust has been detected by seismic studies in the Oceahus 
Procellarum region (ToksiJz et al . , 1972a, b, d) . From the 
orbital x-ray fluorescence experiment (Adler et al. , 1972a, 
b) , from seismic velocities and from analyses of surface 
rocks at or derived from highland sites (Wood et al..-, 1970, 
1971; Turkevich, 1971; Reid et al., 1972), the original 
crust is thought to be . largely feldspathic in composition, 
gabbroic anorthosites with some norite (KREEP) component. 

The chronology of lunar igneous activity starts from 
the formation of the original crust about 4.6 b.y. ago as 
deduced from model ages of lunar rocks and soils (Wasserburg 
and Papanastassiou, 1971; Tatsumoto, 1970; Tatsumoto et al. , 
1972; Silver, 1970) . Although the record of pre-mare vol- 
canism has been obscured by subsequent events, there is 
evidence for activity in the time between the appearance 
of the original crust and the filling of the mare basins. 

The formation ages of KREEP basalts are thought to be about 
4.4 b.y. (Meyer et al^. , 1971; Papanastassiou and Wasserburg, 
1972b; Schonfeld and Meyer, 1972). From these ages and 
Rb-Sr systematics (Papanastassiou and Wasserburg, 1971a) it 
can be concluded that the original, crust differentiated before 
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the formation of KREEP basalts; i.e. within the first few 
hundred million years of lunar history. 

To create a 60-km thick feldspathic crust requires the 
melting and differentiation of at least one half the total 
volume of the Moon (Wood et cd . , 1971). This is an important 
constraint which we will use in this study for the choice of 
the initial lunar temperature profile. 

The igneous activity of the Moon which resulted in the 

filling of mare basins has been well documented by age-dating 

40 39 

of lunar samples, particularly using Rb-Sr and Ar -Ar 
techniques (Papanastassiou and Wasserburg, 1970, 1971a, b, 

1972a; Wasserburg and Papanastassiou, 1971; Turner, 1970, 

1972; Turner et al . , 1972; Husain et ad . , 1972a, b; York 

et ad. , 1972; and others). To date, the well-determined ages 

of mare basalts span a time interval between 3.16 and 3.71 b.y. 

ago (Wasserburg and Papanastassiou, 1971) . Ages for Apollo 14 

rocks are about 3.9 b.y. (Papanastassiou and Wasserburg, 1971b; 

Tatsumoto et al . , 1972; Turner et al . , 1972; Murthy et al. , 

1972) and ages of the crystalline rocks from Apollo 16 are 

greater than 3.8 b.y. (Schaeffer et ad . , 1972; Papanastassiou 

and Wasserburg, 1972b). Along with the 4,1 b.y. age (im- 

placement or reCrystallization) for the sample 15415 (Husain 

et al. , 1972a) and 4.0 b.y. for 12013 (Papanastassiou 

and Wasserburg, 1971b) the above chronology indicates a 

continuum of activity in the Moon from the time of 

formation of the original crust to the time of the 



implacement of mare basalts. The span of igneous activity 
is shown in Figure 1. 

The failure to discover lunar igneous rocks younger 
than about 3.16 b.y. ago is an important constraint on the 
Moon's thermal history (Papanastassiou and Wasserburg, 1971a) . 

If melting occurred since that time, it must have been very 
localized in extent or confined to depths in excess of several 
hundred kilometers. 

2.2 Magnetism of Lunar Rocks 

The lunar crust has an extensive but extremely 

variable magnetization. This has been verified by steady 
magnetic fields measured on the lunar surface (Dyal et al . , 

1970, 1972) and by magnetic anomalies observed by the Apollo 
15 orbital magnetometer (Coleman et al. , 1972). Further, 
many returned lunar samples display a fairly stable remanent 
magnetization (Strangway et al . , 1970; Nagata et al . ,1970 , 

1972; Helsley, 1970; Runcorn et al., 1970; Pearce et al . , 

1972; Collinsori et al . , 1972; Dunn and Fuller, 1972; and others) 
On the basis of these observations, which appear to require 
1000-y fields on the lunar surface during the emplacement of 
rocks spanning a time interval of almost 1 b.y., several 
workers (Runcorn et al., 1970; Pearce et al., 1971, and 
others) have suggested that the fields are internal in origin, 
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and thus indicate a molten conducting core or extensively 
molten pockets that could sustain dynamos between 3.2 and 
4.0 b.y. ago. There is debate over this conclusion/and 
many other explanations of the magnetization, including 
fields associated with the earth, with the solar wind, or 
with shock or thermal effects, have been offered (Hide, 

1972) . 

2.3 Lunar Mascons 

The lunar mascons (Muller and Sjogren, 1968) are 

generally interpreted as arising from higher density mare 

basalts filling the basins excavated in the lower-density, 

anorthositic , initial crust (Sjogren et al. , 1972). The 

time of the mascon formation then can be dated by the ages 

of mare basalts, and it is about 3.3 b.y. for Mare Imbrium 

(Wasserburg and Papanastassiou, 1971). 

Thus the excess masses at mare basins have not been 

totally compensated over a three billion year time span. 

This requires that (1) at the time of the mascon emplacement 

(mare flooding) the lunar lithosphere was thick enough to 

support the loads (Urey, 1968; Solomon and Toksttz, 1972) , and 

(2) the viscosity of the outer layers of the Moon must be 
27 

about 10 poise in order that the stresses have not relaxed 
since mascon formation (Urey and MacDonald, 1971; Baldwin, 
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1971; Arkani-Hamed, 1972). Even in the absence of water 
and volatiles, this high viscosity requires temperatures 
significantly lower than the melting curve in the outer 
several hundred kilometers of the Moon. 

2.4 Electrical Conductivity and Present-Day Temperature 
The electrical conductivity structure and related 
properties of the lunar interior have been deduced from the 
analysis of the Moon's response to transient magnetic field 
events using simultaneous magnetometer measurements on the 
anti-solar surface and in the free-streaming solar wind 
(Dyal and Parkin, 1971; Dyal et al . , 1972), and analysis of 
the frequency-dependent transfer function of sunlit lunar 
surface to the free-stream magnetic field (Sonett et al. , 

1971, 1972). The electrical conductivity profiles were 
compared to laboratory measurements (England et al., 1968; 

Kobayashi and Maruyama, 1971) to show that temperatures 
in the outer 700 to 800 km of the Moon must be less than 
1000°C (Sonett et al. , 1971, 1972; Dyal and Parkin, 1971; 

Dyal et al-.., 1972). When more resistive material is 
assumed for interpreting the electrical conductivity 

(Duba et al., 1972; Anderson and Hanks, 1972; Housley and Morin, 1972) 
temperatures in excess of 1000°C at depths greater than 350 km are 
possible. The suggested temperature limits of Dyal et al. 

(1972) and Duba et al. (1972) are shown in Figure 2. The 



original temperature profile of Sonett et al . (1971) is 
also indicated. 

One property common to all temperature models inferred 
from the electrical conductivity profile is that, at least 
in the outer 600 km, the temperatures are below the solidus 
of anhydrous basalts. This constraint must be satisfied by 
the thermal evolution models. 

2.5 Seismic Velocities and Lunar Tectonism 

The propagation of both compressional and shear waves 
in the lunar mantle requires that to the lowest limit of 
penetration of these waves the mantle must be solid (Toksttz 
et al . , 1972b, d) . From all evidence, moonquakes occur at 
depths as great as 800 km (Latham et al . , 1972a, b) . The 
lunar mantle must be solid at that depth to sustain the 
stresses that cause moonquakes. Nearly constant velocity 
in the upper mantle and the relatively high Q values (Latham 
et al., 1972b) further support the idea that the temperatures 
do not reach the melting point in the upper several hundred 
kilometers of the Moon. 

The low level of tectonic activity on the Moon is also 
consistent with low temperatures and high strength in the 
lunar mantle. Although a number of small moonquakes have 
been recorded by the lunar seismic network, the total seismic 


energy release within the Moon is about eight orders of 
magnitude less than that in the earth (Latham et al. » 1972b). 
The absence of large moonquakes and the small amount of total 
energy release imply the absence of significant dynamic 
activity in the Moon. 

2.6 Lunar Surface Heat Flux 

2 

The Apollo 15 heat flow value of 33 ± 5 ergs/cm -sec 
(Langseth et al . , 1972) places narrow limits on the total 
radioactivity, and to a lesser extent, on the temperatures 
in the Moon. Caution is advisable in adopting this single 
measurement of heat flow as an average value for the Moon, 
especially in light of the variation of KREEP abundances 
c&etected by the orbital y-ray experiment (Metzger et al . , 
1972) . However, the absence of water and of extensive 
tectonic activity argues against strong local variability 
of sub-surface temperature and surface heat flux. The local 
variability of heat flow due to KREEP abundance is probably 
less than a factor of two. Thus we will utilize the measured 
heat flow value as a strong constraint in our thermal models. 



3. MODEL PARAMETERS AND COMPUTATIONAL PROCEDURE 

To calculate theoretically the thermal evolution of 
the Moon, one must specify the initial conditions, the 
sources of heat energy, and the mechanisms of heat transfer. 

The heat sources largely fall into two classes: those that 

affect only the very early history of the Moon and can be 
incorporated in the initial temperature distribution, and 
those that are important throughout lunar history. Heat 
transferred by lattice conduction and radiation can be fully 
accounted for; convective transport must be modeled in an 
approximate fashion. In this section we discuss in detail 
the necessary physical parameters and the particulars of our 
computational scheme. 

3.1 Initial Temperature 

The initial temperature distribution in the Moon is a 
function of the process of lunar formation and of the immediate 
environment of the Moon during and shortly after its origin. 

The evidence favoring an early episode of extensive near- 
surface melting in the Moon has been outlined above . Energy 
sources which might have heated the Moon to solidus or near- 
solidus temperatures include gravitational energy liberated 
during accretion, tidal dissipation, solar wind flux, short- 
lived radioactivity and adiabatic compression. 



The latter two heat sources are likely of small con- 
sequence. The energy released by adiabatic compression in 
the Moon does not raise the temperature by more than a few 
tens of degrees. The effect of short-lived radioactive 
isotopes depends on the difference in time between nucleo- 

genesis and the accretion of the Moon. The work of Schramm 

26 

et al . (1970) has demonstrated that A1 , previously thought 
to be perhaps the most likely radioactive source of substantial 
heat in the early history of the solar system (Fish et aX. ,- 
1960) ,. was not an important heat source in meteorites or in 
the Moon. 

Both tidal dissipation and solar-wind flux can, under 
favorable circumstances, heat most or all of the Moon to 
melting temperatures. The effect on lunar temperature of 
tidal dissipation is small (Kaula, 1966) except when the 
Earth-Moon distance is only a few earth radii. This might 
have occurred during a short time interval after a capture 
event or after the Moon formed in close proximity to the earth 
(Singer, 1970). How this heat would be distributed with depth 
is uncertain, though we might qualitatively expect that most 
of the anelastic losses occur in those portions of the Moon 
that are at the highest temperatures before the tidal heating 
episode began. Heating of the Moon by the electric currents 
produced by a uni-polar generator driven by the Solar wind 
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might be sufficient to melt the Moon if the Sun passed 
through a T-Tauri stage (Sonett and Colburn, 1970) , otherwise 
the effect would be small. The likelihood that either of 
these energy sources was important is difficult to assess. 
Neither is explicitly included in our models, but recognition 
is made of the possibility that virtually the entire Moon 
could have melted early in its history. 

Since it is reasonable to assume that the Moon formed 
by some accretion process, the most promising initial heat 
source for providing energy sufficient to melt the near- 
surface regions of the Moon is gravitational energy of 
accretion. Suppose that at some time t during accretion 
the proto-moon has radius r, and that during a time increment 
dt an outer shell of thickness dr is added to the lunar 
mass. Then the temperature T within that shell may be 
estimated by equating the additional gravitational potential 
energy to radiated energy plus heat: 

pGM|r)_ dr = ea (T 4 - T fa 4 )dt + pC p (T - T b )dr (1) 

where p is the density of the accreting particles, G is 
the gravitational constant, M(r) is the mass of the body 
within radius r, a is the Stefan-Boltzman constant, e is 
the emissivity (here. taken equal to 1), and C p is the 



specific heat. is the base temperature of the accreting 

particles, or the equilibrium temperature of the cloud of 
material within which the Moon formed. It will reflect 
whether the Moon formed as an aggregation of cold particles 
or condensed from a cooling gas; any effects of short-lived 
radioactivity can be included. T(r) can be calculated from 
Equation (1) by specifying p, C p , and the accretion 

rate dr/dt. An ad hoc but physically reasonable model for 
accretion was suggested by Hanks and Anderson (1969) , and 
it has been adopted in this paper (Figure 3a) . The accretion 
rate is expected to be slow when the Moon is just beginning 
to form, to increase as the accreted mass increases and then 
to taper off to zero as the material is depleted and the 
-final lunar radius is achieved. These ideas are adequately 
represented by an accretion rate 

dr/dt = ct sin y t. (2) 

The constants c and y are determined from the Moon's 
final radius and the total duration of accretion. 

If the accretion time is considered arbitrary, then it 
may be adjusted so as to give an initial temperature profile 
in excess of the solidus in the outer portions of the Moon. 
Using Equations (1) and (2), it can be shown that the near- 



surface temperatures never exceed the basalt melting curve 
for accretion times greater than about 1000 years. This 
conclusion was also reached by Mizutani et al. (1972) from 
a more detailed analysis of the accretion process. One 
can estimate the duration of accretion by matching this model 
to the temperatures required to cause melting to specified 
depth. To differentiate a 60 -km crust composed largely of 
feldspar probably requires melting to a depth of at least 
500 km; i.e. , involving at least half of the Moon's volume. 

This can be achieved if the total accretion time is about 
100 years. The duration of accretion may be made arbitrarily 
long, however, by appropriately lowering the assumed emissivity 
of the lunar surface. The temperature profile resulting from 
Equations (1) and (2) for a 100-year accretion time is given 
in Figure 3b for two separate base temperatures; these pro- 
files are adopted as starting temperatures in the calculations 
below. 

3 . 2 Heat Sources 

The most important heat source in the thermal history 
calculations is the heat generated by the long half-life 
radioactive isotopes , Th^^ , and . The abundances 

of these isotopes measured in many lunar samples and the 
orbital y-ray spectrometer results (Metzger et al. , 1972) will 



permit surface concentrations to be estimated over much of 
the Moon. While we can certainly conclude that the abundances 
at the Moon's surface must be considerably higher than at 
depth (Papanastassiou and Wasserburg, 1971a; and others) , 
the bulk radioactivity of the Moon and the distribution of 
the radioactive isotopes with space and with time, are difficult 
to constrain because of the complexities introduced by 
melting and crystal-liquid fractionation processes. To 
guide our estimates we consider the radioactivity measured 
in lunar samples and meteorites and the determination of heat 
flow at the Apollo 15 site. 

In Figure 4 , total K and U abundances are plotted 
for Apollo samples and some terrestrial basalts and meteorites . 
This figure is an updated version of one given by Hays 
(1972). For all lunar rocks (including basalts, breccias 
and soils) a constant K/U ratio of about 2000 is a good 
approximation although there are small, systematic departures 
from this average and the absolute abundances vary. This 
value compares with K/U = 10,000 for the earth and 80,000 
for chondrites. The low K/U ratio for the Moon is linked 
to a general depletion of volatile elements (LSPET , 1969; 

1972). Clearly neither terrestrial nor chondri tic abundances 
should be used for lunar calculations. Chondrites can also 
be rejected on the basis of Rb abundances as a model for the 
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primitive lunar material (Papanastassiou and Wasserburg, 
1971a). The Th/U ratio for lunar rocks is about 3.6 to 
4.0, consistent with terrestrial and meteorite values; 

As a possible handle on the average abundances of radio- 
active isotopes in the Moon we consider the achondrites, 
especially howardites and eucrites. A genetic relationship 
between achondrites and the Moon has been suggested by sev- 
eral investigators (Duke and Silver, 1967; and others) ; On 
the average, howardites have a U abundance of about .035 
ppm while for (brecciated) eucrites, which are enriched in 
refractory elements, the value is about .130 ppm (Mason, 

1971) . It is of interest that the K/U ratio of eucrites 
is close to that for lunar rocks, while that of howardites 
(based on fewer samples) is somewhat higher. 

In this paper we adopt the scaling of radioactive 
elements K/U = 2000 and Th/U = 4 throughout the present- 
day Moon. The average radioactivity of the Moon is then a 
single- valued function of the bulk uranium content. We con- 
sider models with U concentrations appropriate to howardites 

and eucrites and also concentrations set so as to match the 

2 

Apollo 15 heat flow value of 33 erg/cm -sec (Langseth et al . , 

1972) . The relevant decay constants, decay energies and iso- 
topic abundances are taken from Clark (1966) . 

We should note that we ignore tidal dissipation as a 
possible long-term heat source of some significance. This 
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is done because of both the difficulty in quantitatively 
modeling the phenomenon and the uncertainty in past lunar 
orbit characteristics. As a result, we may be underestima- 
ting the heat generation in the early history of the lunar 
interior . 


3.3 Thermal Conductivity and Other Parameters 

In the absence of convective flow, transport of heat 
is governed by an effective thermal conductivity equal to 
the sum of lattice conduction and radiation terms. The 
effective thermal conductivity for some lunar and terrestrial 
materials as well as theoretical estimates (MacDonald, 1963b; 
Schatz and Simmons, 1972) are shown in Figure 5. Clearly, 
the values for lunar mare basalts should not be used since 
these rocks do not represent the bulk of lunar composition. 
(Ringwood and Essene, 1970; Gast, 1972) . In this study we 
adopt the curve proposed by Schatz and Simmons (1972) for a 
hypothetical dunite in which radiative transfer is somewhat 
suppressed by grain-boundary extinction. At low temperatures, 
the conductivity is primarily due to lattice conduction and 
it decreases with increasing temperatures. At higher tem- 
peratures, the radiative term, roughly linear with tempera- 
ture, is dominant. Previous studies of thermal evolution of 
the Moon have generally used MacDonald's formulation and 
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assumed that the opacity is constant with temperature. This 
gives the radiative contribution a cubic dependence on tem- 
perature. The pressure effects are probably small at lunar 
pressures and are neglected. In the calculations an analytic 
expression has been used for the curve marked "S" and con- 
ductivity is adjusted at each step for appropriate temperature. 

The melting curve, surface temperature, specific heat, 
heat of fusion, and density must be specified for calculations. 
The melting curve used in this paper is the solidus from the 
melting range of mare basalts (Ringwood and Essene, 1970). 

We used only one curve to specify both partial and total 

3 

melting. A constant density of 3.34 g/cm and a constant 
surface temperature of -20°C (Langseth et al., 1972) have 
been used. A summary Of the parameters common to all models 
is given in Table I. 

3.4 Computational Technique 

The thermal calculations are carried out for a spheri- 
cally symmetric moon (parameters varying with radius only 
at a given time) , taking into account melting, simulated 
convection and differentiation. In this section these 
steps are described briefly. 

Temperature evolution models are calculated using the 
finite-difference solution of the heat conduction equation 
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where C p is the specific heat, p is the density, T is 
the temperature, r is the radius, K is the thermal con- 
ductivity and H(r,t) is the heat-source term, which in 
general depends on radius and on time. K is taken to be 
temperature-dependent as discussed above. 

The finite-difference analog of Equation (3) which con- 
serves heat flux is 
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where 


p = Ar 
r = np 
t = mAt 


a = At/C p p 
T™ = T(np,mAt) 


Since the scheme is explicit, a stability condition relating 
the time increment and the grid spacing must be fulfilled. 
Questions of stability, convergence and accuracy of Equation 



(4) are discussed fully by Toksflz et aly (1972c). The 
stability condition which is developed in the appendix of 
that paper is 


At < 


C p pAr 


n 


(n ■+ 1/2) 



where K is the maximum value of the conductivity at a 
max . 

given time step. The time increment is computed at each 
time step using Equation (5) with a factor of 4 instead of 
2 in the denominator. Generally a grid spacing of 20 km is 
used . ' 

To avoid a temperature singularity at the center of the 
Moon (r =0), the solution was started at r = 20 km. 

Heat flux was taken as zero at the grid point corresponding 
to r - 20 km by setting T^ = T^ where T^ is the tem- 
perature at r^ = 20 km and T^ is the temperature at 
r £ = r^ + r. Temperature was held constant (-20°C) at the 
surface of the Moon. Surface heat flux was computed by 
adding the heat flow at the boundary between the upper two 
finite-difference cells to the steady-state heat flow due 
to heat production in the uppermost half-cell; the former 
quantity was corrected for the slightly greater area of 
the lunar surface. 
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Simulation of Melting and Convection . The effects of melting 
and convection on temperature were modeled using a technique 
described by Reynolds et al. (1966). In this technique, the 
temperature profile at some time mAt is compared with a 
given melting-point curve. If the temperature is above 

the melting point, 'T , the temperature increment above the 
melting point (T™ - T r ) is converted to its heat equivalent, 
AH, by dividing by the specific heat and the density. If 
AH is not equal to or greater than a specified heat of 
fusion, the material is taken as partially molten, the tem- 
perature at point nAr is held at T n and another iteration 
in time performed. The material at point nAr remains par- 
tially molten at temperature T^ until the sum of the AH's 
from mAt to some (m + p)At equals the heat of fusion. 

At this point in time the material is completely molten and 
temperature is allowed to increase above T .• Provision is 
made for the phase change to run either way so that molten 
material may solidify by releasing heat equal to this heat 
of fusion. 

Simulation of convection of molten material follows 
the technique outlined above except that once the material 
has become completely molten, temperature is held at the 
melting point. Any increase in temperature which would 
raise the temperature above the melting point is converted 



into its heat equivalent and transferred upward to the next 
grid point as an equivalent temperature increment due to 
convection. The conversion to equivalent heat and back 
to temperature is necessary since the volume element in- 
creases with radius. 

The main assumption underlying this technique is that 
convection occurs on complete melting and limits the temper- 
ature in the region of complete melting to the melting tem- 
perature. This scheme accounts for the transfer of heat 
which would have taken place under convection without actu- 
ally using the mass transfer terms in the equation. Con- 
vection in a partial melt can be approximated by adopting 
an appropriately lower heat of fusion. Convection by solid- 
state creep can be roughly approximated by redefining the 
"solidus" (see below). 

Differentiation . In the regions where melting occurs , the 
magma would be enriched with U, Th and K and eventually 
would transfer these heat sources toward the surface. To 
account for this, at discrete time steps all heat sources 
to the depth of deepest melting are differentiated upwards. 
Where no melting has occurred, no differentiation is carried 
out and initial radioactive abundances are maintained. 

The concentration of radioactive heat sources after 



each differentiation is assumed to decrease exponentially 
with depth (Lachenbruch, 1968) : 


H (r) 


A exp 
o 



(R- 

h 




(6) 


where A is the surface abundance and R is the lunar 

o ••• 

radius. The skin depth h is presumed to decrease with 
time. The skin depth of the last differentiation is set 
so that the present-day radioactivity at the lunar surface 
equals that of the average of soil samples from mare sites 
about 1 ppm (Silver, 1970; 1972; and Figure 2) , with the 
exception that the skin depth is never taken less than the 
finite-difference grid size. Differentiation is carried 
out in 1 to 4 steps, starting with h in the range 100 to 
200 km and ending with h = 20 to 60 km. The time of last 
differentiation is taken to lie in the range 3.0 to 3.7 
b.y. ago. Thermal models calculated by these prodedures 
are discussed in the next section. 
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4. TEMPERATURE MODELS 

In this section we present specific models for the thermal 
evolution of the Moon. An assortment of initial conditions 
and bulk concentrations of radioactive isotopes are employed 
both to demonstrate the effect such parameters have on thermal 
history and to assess whether the constraints discussed in 
Secion 2 can limit models of lunar formation or constitution. 

A quick overview of all models presented below is given in 
Ta±>le II . 

We begin by considering the evolution in a uniform, ini- 
tially cold (T = 0°C everywhere) Moon. While such a model 
is unrealistic, it is a standard output of those presenting 
thermal history calculations and serves to verify the compu- 
tational scheme. In Figure 6 is shown the present-day tem- 
perature profile in an initially cold Moon as a function of 
the present-day concentration or uranium. For U 5 37 ppb 
(parts per billion), the temperature is below the basalt 
solidus throughout the Moon and has been so since the Moon's 
formation. Clearly models such as those in Figure 6 that 
begin with an eyerywhere-cold moon cannot generate enough 
heat from long-lived radioactivity alone to provide the ne- 
cessary melting and differentiation in the Moon's upper few 
hundred kilometers. ' 
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As a model with an antithetical initial state, we show 
the thermal evolution in an initially molten Moon in Figure 
7. As mentioned in Section 3, complete melting in the Moon 
shortly after its formation could have been produced by 
close approach of the Moon to the Earth or by immersing the 
Moon in the solar-wind flux of a Sun passing through a 
T-Tauri phase. Upon complete melting in the model, all the 
radioactive elements are concentrated near the lunar surface. 
In such a model, iron or dense iron compounds (such as FeS) 
would differentiate to the Moon's center to form a core. 
Because of the depletion of heat sources in the interior 
and because of the efficient heat transfer in molten re- 
gions, the Moon cools relatively rapidly and is everywhere 
solid after 2 billion years. At present the Moon is cooling 
throughout; the average concentration °f uranium is 60 ppb 
(corresponding, perhaps, to a mixture of howardite and eucrite 

2 

compositions) , and the surface heat flow is 29 erg/ cm -sec, 
within the range of uncertainty of the Apollo 15 measure- 
ment. A similar model with less radioactivity was given in 
ToksGz et al . (1972c). 

Though early melting throughout the Moon can probably 
be produced only by a catastrophic event, the normal pro- 
cess of accretion can yield melting of the near-surface 
regions sufficient to segregate an early crust and to 



78 


differentiate the radioactive elements. We now consider a 
series of thermal history models with one of the initial 
temperature profiles from Figure 3. The primary variable 
of interest will be the concentration of radioactive heat 
sources . 

A thermal evolution model with average present-day 
uranium concentration equal to that of howardites (35 ppb) 
is shown in Figure 8. The initial temperature profile is 
derived from Equations (1) and (2) with a base temperature 
of 500°C and a total accretion time of 100 years. Melting 
initially extends to 600 km depth, and the interior is rela- 
tively cool. In this model, the outer portions of the Moon 
cool monotonically with time while the interior progressively 
warms. Between 0.5 and 1.5 b.y. after lunar origin, a solid 
lithosphere on the Moon grows from 150 to 310 km in thick- 
ness. Partial melting occurs at some depth over the entire 
lunar history, though upwards differentiation of radioactivity 
was arbitrarily stopped in the model 3.7 b.y. ago. The inner 
950 km is partially or completely molten at present, with 
heat coming from undifferentiated material below 900 km 

depth. The average heat flow at the Moon's surface in the 

2 

model is currently 18 erg/cm -sec. 

In Figure 9 is given a thermal evolution model with 
the same initial temperature profile as in Figure 8 but 


with an average present-day concentration of uranium equal 
to that of eucrites (130 ppb) . In this model, the entire 
Moon below 50 km depth becomes completely or partially 
molten within 0.5 b.y. after lunar formation. The lunar 
lithosphere grows very slowly; at 1.5 b.y. it is only 70 km 
thick. Such a value is probably much too small if the 
mascon-generated stresses were supported in the litho- 
sphere. At the present time in the model, the Moon is 
partially or completely molten below 200 km depth, except 

for a small solid core, and the surface heat flow is 55 erg/ 

2 

cm -sec. Models similar to that of Figure 9 have been pro- 
posed elsewhere (Hanks and Anderson, 1972) . 

We feel that the thin lithosphere in the first 3 b.y. of lunar history, 

the high heat- flow, and the difficulty of explaining moon- 

quakes 700 to 800 km deep make such a model an unlikely 

candidate. While eucrites might be linked to lunar rocks, 

the average radioactivity of the Moon cannot be as high as 

that of eucrites. 

Since the present-day heat flow scales approximately 

as average uranium concentration in thermal history models 

in which other parameters are held fixed (Solomon and Toksflz, 

1972), we can determine the radioactivite element content 

necessary to match the single measurement of heat flow at 

2 

the Moon's surface, 33 ± 5 erg/cm -sec (Langseth et al . , 



1972). We should caution that the y-ray spectrometer re- 
sults (Metzger et ad., 1972) indicate that the region sur- 
rounding the Apollo 15 site is significantly more radio- 
active than most of the lunar surface,- so the average heat 
flow on the Moon may be less than the Apollo 15 value. 

A model giving a present-day surface heat flow of 
33 erg/cm -sec is shown in Figure 10. The initial temper- 
ature profile is from the accretion model of Figure 3; a base 
temperature of 0°C and a total accretion time of 100 years 
were used. The last differentiation of radioactive heat 
sources was taken to occur 3.7 b.y. ago. The present-day 
uranium concentration averages to 70 ppb. The time sequence 
of melting in the model is shown. in Figure 11. 

It may be seen that the model satisfies the major con- 
straints upon lunar evolution. The zone of melting deepens 
with time. The lithosphere during the first 2 b.y. thickens 
at the rate of about 160 km/b.y. In particular, the shallow- 
est melting progresses from 155 to 245 km depth during the 
period of mare filling, in agreement with the depth of origin 
of mare basalts (Ringwood and Essene, 1970) and with the 
need for a reasonably thick lithosphere to sustain the 
stresses associated with mascon gravity anomalies. Between 
1 and 2 b.y. in the model, a convecting core develops due 
to the undifferentiated and relatively radioactive primordial 


material near the Moon’s center. The convecting core grows 
to a radius of 1200 km at 3 b.y. and then slowly shrinks as 
the concentration of heat sources decays with time. The 
Moon is molten below 580 km depth at present in the model. 
This is in apparent conflict with the occurrence of moon- 
quakes at 700 to 800 km depth, unless such events are not 
representative of the average Moon, and implies that a 
Moon with somewhat less radioactivity, and somewhat lower 
heat flow, is more appropriate . The feature of a molten 
convecting core is absent in models (such as those of 
Figures 7 and 9) in which the Moon differentiates completely 
and concentrates all of the radioactive heat sources toward 
the surface . 

Models similar to that of Figure 10 except starting 
at a higher base temperature and . involving more recent 
segregation of radioactive heat sources from the partially 
molten regions of the Moon can match the Apollo 15 heat 
flow value with a somewhat lower average uranium concentra- 
tion. With a range of initial temperature profiles from 

accretion models, a present-day surface heat flow of 33 ± 

2 

5 erg/cm -sec implies an average uranium concentration in 
an originally homogeneous Moon of 65 ± 15 ppb (Solomon and 
Toksflz, 1972) . 



The hypothesis of a stable magnetic dynamo in a 
liquid, iron-rich lunar core from 3 to 4 b.y. ago would re- 
quire the deep interior of the Moon to begin at temper- 
atures at least in excess of the Fe-FeS eutectic or to 
have formed at lower temperatures with substantial 
radioactivity at depth (cf. Brett, 1972). The thermal 
history models of Figures 7 and 10 are consistent with this 
idea. In Figure 10, temperatures .5 b.y. after lunar 
origin exceed the Fe-FeS eutectic (Brett and Bell, 1969) . 
throughout most of the Moon. Thus core formation by 
melting, sinking and aggregation of iron-rich droplets is 
achieved prior to the formation of the oldest lunar rocks 
possessing remanent magnetization. 

Except for models which are initially molten and 
segregate all radioactivity upwards toward the surface, thermal 
history models beginning with uniform distribution of radioactivity 
and satisfying the Apollo 15 heat flow value are presently 
partially molten or at temperatures within 100°C of the 
solidus at depths below about 500 km (Solomon and ToksGz, 

1972). Thus for such models the heat flow and deep moonquakes 
are difficult to reconcile. One means of avoiding the 
dilemma is to make an additional postulate: the concentration 

of radioactivity in the Moon decreases with depth as a 
primary feature associated with chemically inhomogeneous 
accretion, as suggested by several investigators 
(Gast, 1972; Green et al., 1972; and others) . 



Shown in Figures 12 and 13 is a thermal history model fol- 
lowing closely one proposed by McConnell and Gast (1972) . 

The radioactivity is stratified: if no differentiation 

had been allowed, the present-day abundance of uranium 
would be 120 ppb in the uppermost 200 km, 50 ppb in the 
next 250 km, and 10 ppb in the remainder of the Moon. . 

The average present-day uranium concentration is thus 53 
ppb. The ratios K/U and Th/U are held constant. The 
initial temperature profile is from the accretion model of 
Hanks and Anderson (1969) with a base temperature of 0°C and 
an accretion time of 1000 years. The initial temperature 
barely exceeds the solidus between 40 and 140 km depth; a 
shorter accretion time would have raised the initial tem- 
perature. Differentiation of radioactive elements is 
accomplished in four discrete steps between 0 and .9 b.y. 

The temperature evolution in the upper few hundred kilo- 
meters of the Moon is similar to models of comparable bulk 

radioactivity discussed earlier. The present-day heat flow 
2 

is 32 erg/cm sec and the Moon is everywhere at least 200°C 
below the solidus . 

The difference in the behavior of the model of Figure 
12 from those discussed earlier is that the center of the 
Moon in this model is tacitly assumed to have been relatively 
depleted in radioactive elements at cold temperatures. Thus 
the lunar center heats slowly and the Moon is everywhere 
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solid at present. Whether a cold accretion can accompany 
primary chemical layering is an important question. Such 
an origin is incompatible with a molten iron-rich core at 
an early stage in lunar history. If questions such as these 
can be resolved and if the Apollo 15 heat flow value holds 
up, then the temperatures in the deep lunar interior, if 
discernible, offer the only possible means by which thermal 
history calculations can distinguish primary accretiohal 
zoning from extensive in situ differentiation. 

In all of the models above, heat transfer by convective 
mass transport was ignored at sub-solidus temperatures. This 
may be quite unreasonable, particularly for the deep lunar 
interior. Several authors (Runcorn, 1962; Kopal, 1962; 
Turcotte and Oxburgh, 1969b; Tozer, 1972) have argued that 
convection by solid-state creep should be at least as 
efficient a heat-transfer mechanism as lattice conduction 
and radiation at temperatures considerably below the solidus. 
Inclusion of the effect of solid-state creep should produce 
present-day temperature profiles that are cooler, than those 
given above. Tozer (1972) , in fact, has stated that the 
relatively cool lunar temperature profiles proposed by Sonett 
et al . (1971) and others interpreting the electrical conduc- 
tivity distribution are a natural consequence of a convecting 


interior . 
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While we do not wish to attempt to solve the complete 
problem of thermal evolution in a realistically convecting 
Moon, we can make some simplifying assumptions that allow 
us to make a rough estimate of the effect of solid-state 
creep on our thermal history calculations. Convective 
solutions to the present day temperature profile in the 
Moon generally show a stable shell of a few hundred kilo- 
meters overlying a convecting interior of nearly constant 
viscosity (Turcotte and Oxburgh, 1969b; Tozer, 1972), a 
consequence of the low pressures and small adiabatic temper- 
ature gradient in the Moon. We therefore postulate in our 
thermal evolution models that solid-state creep will be an 
efficient heat transfer mechanism when the temperature i n 
the bulk of the lunar interior is such that the viscosity 
exceeds a critical value. This temperature is essentially 
the 'stabilization temperature' of Tozer (1970, 1972). Con- 
vection is simulated, as before, as a 100-percent-efficient 
process: all heat in excess of that required to maintain the 

viscosity of the material at the critical value is transferred 
upward in the model. We assume a rheology proposed by 
Turcotte and Oxburgh (1969a) for diffusion creep in the earth's 
mantle; whether or not this is compatible with our assumed 
solidus is debatable. Thus we ignore the probable complica- 
tion of non-Newtonian viscosity (Weertman, 1970) , an uncertain 



increase in viscosity associated with the depletion of 

volatiles in the Moon relative to the earth (Orowan, 1965), 

and the apparent conflict between geophysical constraints 

on the rheology of the earth 1 s mantle and presently available 

laboratory measurements of high-temperature creep in rocks 

(Goetze and Brace, 1972). For the critical value of vis- 

24 

cosity we chose 10 poise. This value is an estimate, 
obtained from marginal stability theory in a uniform 
gravitating sphere, of the maximum value permissible if the 
Moon is currently convecting (Turcotte and Oxburgh, 1969b). 

We found that using lower values for the critical viscosity 
did not alter significantly the thermal evolution from models 
similar in all respects but neglecting solid-state creep? 
this is because the temperatures at which such viscosities 
are reached, according to the viscosity-temperature relation- 
ship of Turcotte and Oxburgh (1969a), are near the basalt 

2i 

solidus. Tozer (1972) has proposed that a viscosity near 10 
poise and less creep resistant rheologies are more appropriate. 

A thermal evolution model following such a scheme is 
shown in Figure 14. This model begins with an accretion- 
generated initial temperature curve from Figure 3 , and is 
treated as all earlier models for the first 2 billion years. 

At that time, almost all of the Moon below 250 km depth is 
at a temperature high enough to give a viscosity less than 



24 

10 poise. Although the viscosity was less than that 
value over a more limited depth range at times less than 
2 b.y., we argue that solid-state creep was probably un- 
important at these early stages by appealing qualitatively 
to the very strong dependence of the Rayleigh number, a 
criterion for hydrodynamic instability, on the characteristic 
length (depth interval) of the convective cell (e.g. Turcotte 
and Oxburgh, 1969b). After 2 b.y. in the model, the temper- 
ature in the convecting interior is stabilized at that 

24 

temperature necessary to maintain a viscosity of 10 poise; 

the gradient is slightly superadiabatic . All excess heat 

is transferred to the non-convec ting outer shell, which is 

slowly cooling and thickening with time. The present-day 

average uranium concentration is 23 ppb and the heat flow 

2 

at the lunar surface in this model is 13 erg/cm sec . For 
a higher uranium concentration, the heat flow will be 
higher and the only change in the present-day temperature 
profile from that given in Figure 14 will be a slightly 
thinner lithosphere . 

All of the thermal models discussed above and summarized 
in Table II demonstrate the importance of the initial conditions 
and composition for the evolution and present-day properties 
of the Moon. With additional data, it is probable that the 
early history of the Moon can be reasonably well understood 
from computations such as these . 
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5. CONCLUSIONS 

A wide assortment of data now constrain models of the 
evolution of the Moon's interior. We have quantitatively 
examined some of these constraints and the suite of lunar 
thermal models consistent with them. From these we can 
conclude : 

(1) Models of the thermal evolution of the Moon that fit 
the chronology of surface igneous activity, the need 
for a cool, rigid lithosphere since the formation of 
mascons, and the surface concentrations of radioactive 
elements all require extensive differentiation and 
upwards concentration of radioactive heat sources early 
in lunar history. This differentiation may be either 

a primary accretional feature or a product of nearly 
Mbon-wide melting. 

(2) A model of the Moon which is initially hot, and exten- 
sively or totally molten can satisfy most of the con- 
straints described in Section 2. A totally molten Moon 
cools faster than other models due to complete and early 
differentiation . 

(3) The early history Of the Moon is characterized by strong 
magmatic and tectonic activity. Convection must have 
played an important role in lunar differentiation and 
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dynamics. Patterns of convective flow will introduce 
lateral heterogeneities which cannot be modeled with 
a radially symmetric Moon. 

2 

(4) If the Apollo 15 heat flow value of 33 ±5 erg/cm -sec 
(Langseth et al . , 1972) is representative of the Moon, 
then the average concentration of uranium in the Moon 
is currently 65 ± 15 ppb, assuming the Th/U and K/U 
ratios of lunar surface rocks hold throughout the 
interior. This uranium concentration is intermediate 
between the average values for howardites and 
eucrites, although closer to howardites. Such an 
achondritic composition would satisfy the bulk of the 
geochemical data as well as the mass and moment of inertia 

. constraints for the Moon. 

(5) Thermal history calculations are only marginally useful 
for distinguishing between primary accretional zoning 

o>f radioactive heat sources and an early melting episode 
including extensive differentiation. The major difference 
between the two sets of models is the temperature of the 
innermost 1000 km, a region of the Moon for which inde- 
pendent information from seismic velocity or electrical 
conductivity studies is currently lacking. 

(6) Better understanding of the thermal evolution of the 



Moon requires further constraints on heat sources , 
improved computational techniques for convection 
models, and additional data limiting present-day 
internal temperatures. 
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Table I 


Parameters Used in 

Radius 

Density 

Thermal Conductivity 
Heat of fusion 
Specific heat 
Surface temperature 
K/U ratio 
Th/U ratio 
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Thermal Models 

1740 km 
3.34 g/crn* 

Curve "S" in Figure 5 
400 joules/g 
1.2 joules/g°C 
-20°C 
2000 
4 



A Summary of Parameters and Results for Thermal History Models 
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FIGURE CAPTIONS 

FIGURE 1. A summary of igneous activity at the lunar sur- 

39 40 

face from Rb-Sr and Ar -Ar ages of re- 
turned lunar samples. Samples from different 
missions are designated by "A" for Apollo 
(A-il is Apollo 11) and L-16 for Luna-16. Special 
rocks are designated by numbers. Sources for 
this summary are given in full in the text. 

FIGURE 2. Some estimates of present-day temperatures in the 
Moon using the lunar electrical conductivity pro- 
file. Included are the temperature curve of 
Sonett et al. (1971) , an interpretation of the 
Sonett et al^ (1971) conductivity distribution 

using the temperature-conductivity relation of - 

3+ 

an olivine with little or no Fe (Duba et al. , 
1972), and the estimates by Dyal et al. (1972) 
of bounds oh temperature from their most recent 
conductivity models. 

FIGURE 3. a) Two proposed models for the increase in lunar 
radius r with time t during accretion (from 
Hanks and Anderson, 1969) . x is the total ac- 
cretion time. 



b) Initial temperature profiles for a Moon ac- 
creting with a time-dependent radius growth and 
with base temperatures of 0° and 500°C. Total 
accretion time is 100 years. The solidus of 
mare basalt (Ringwood and Essene, 1970) is 
assumed to be an upper limit for possible ini- 
tial temperatures. . 

FIGURE 4. Potassium and uranium concentrations in selected 
meteorites and lunar and terrestrial rocks, up- 
dated from Hays (1972). Points for average chon- 
drites, eucrites and howardites are from Mason 
(1971). Values for Apollo 15 and 16 rocks are 
from LS PET (1972) and LSPET (1973), respectively. 


FIGURE 5. Thermal conductivity of selected materials. 

Sources of measured data (solid lines) are Cremers 

(1971) for Apollo 11 fines, sample density 1.64 

g/cm^; Horai et al . (1970) for Apollo 11 basalt; 

Murase and McBirney (1970) for synthetic Apollo 

11 basalt; Schatz and Sinunons (1972) for sintered, 

polycrystalline forsterite and for single-crystal 

olivine (Fo 0 ,. Fa..) and enstatite. Theoretical 
86 14 

(dashed) curves include those of MacDonald (1963b) 



for two different values of the mean extinction 
coefficient or opacity e, and one proposed by 
Schatz and Simmons (1972) for polycrystalline 
olivine qf approximate composition Fo^q Fa l0‘ 

The Schatz and Simmons curve, labeled S above, 
was adopted in this paper. 

FIGURE 6. Present-day temperature profiles in an initially- 
cold Moon (0°C at all depths) as a function of 
present-day concentration of uranium. The solidus 
of anhydrous mare basalt (Ringwood and Essene, 
1970) is also shown. 

FIGURE 7. A model of the evolution of the lunar temperature 
profile aS a function of time for an initially 
molten Moon. See the text for other parameters. 
Time, in billions of years since lunar origin, 
is indicated by the number adjacent to each pro- 
file. On this and later figures, the Moon is 
partially or completely molten at those depths 
where the temperature profile lies along the soli- 
dus of mare basalt (Ringwood and Essene, 1970); 
the depth range for complete melting is delimited 
by the small arrows above the solidus. 



FIGURE 


FIGURE 


FIGURE 


FIGURE 


8. Thermal evolution in a M6on accreted in 100 years 
at a base temperature of 500°C with average 
present-day uranium concentration equal to that, 
of howardites. All symbols are described in 

the caption to Figure 7. 

9. Thermal evolution in a Moon similar to that of 
Figure 8 except with average present-day uranium 
concentration equal to that of eucrites. All 
symbols are described in the caption to Figure 7 . 

10. Thermal evolution in a Moon accreted in 100 years 
at a base temperature of 0°C with average present 
day uranium concentration chosen to match the 
Apollo 15 heat flow value (Langseth et al ♦ , 1972) 
Shown also are the Fe-FeS eutectic temperature 
(Brett and Bell, 1969) and the phase boundary 
(dotted line) between feldspathic and spinel 
pyroxenite in the model lunar mantle composition 
of Ringwood and Essene (1970). Other symbols 
are described in the caption to Figure 7. 

11. The extent of melting in the thermal model of 
Figure 10. Heavy stippling denotes complete 
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fusion; lighter stippling indicates partial 
fusion. 

FIGURE 12. Thermal evolution in a Moon derived from inhomo- 
geneous accretion. The initial temperature pro- 
file follows from a presumed cold (0°C) accretion 
over a 1000 year time interval. Primary zoning 
of radioactivity is assumed (see text for details). 
Other parameters and symbols are explained in 
the caption to Figure 7. 

FIGURE 13. The extent of melting as a function of time in 

the thermal model of Figure 12. Shading follows 
the notation of Figure 11. 

FIGURE 14. The effect of convection by solid-state creep 

on the later stages of thermal evolution in the 
Moon (see text). From 0 to 2 b.y. , this 
model follows a conventional thermal history 
route beginning with a 100-year accretion at 
0°C. 
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ABSTRACT 

The composition/ structure and evolution of the Moon's 
interior are narrowly constrained by a large assortment of 
physical and chemical data. Models of the thermal evolution 
of the Moon that fit the chronology of igneous activity on 
the lunar surface, the stress history of the lunar lithosphere 
implied by the presence of mascons, and the surface concen- 
trations of radioactive elements, involve extensive differ- 
entiation early in lunar history. This differentiation may 
be the result of rapid accretion and large-scale melting or 
of primary chemical layering during accretion; differences in 
present-day temperatures for these two possibilities are 
significant only in the inner 1000 km of the Moon and may not 
be resolvable. If the Apollo 15 heat flow result is repre- 
sentative of the Moon, the average uranium concentration in 
the Moon is 50 to 80 ppb. 

Density models for the Moon, including the effects of 
temperature and pressure, can be made to satisfy the mass and 
moment of inertia of the Moon and the presence of a low 
density crust inferred from seismic refraction studies only 
if the lunar mantle is chemically or mineralogically inhomo- 
geneous. The upper mantle must exceed the density of the 
lower mantle at similar conditions by at least 5 percent. 

The average mantle density is that of a pyroxenite or olivine 
pyroxenite, though the density of the upper mantle may exceed 

3 

3.5 g/cm . The density of the lower mantle is less than that 



of the combined crust and upper mantle at similar temperature 
and pressure, thus reinforcing arguments for early Moon-wide 
differentiation of both major and minor elements. The 
suggested density inversion is gravitationally unstable and 
implies stresses in the mantle 3 to 5 times those associated 
with the lunar gravitational field, a difficulty that can 
be explained or avoided by (i) adopting lower values for the 
moment of inertia and/or crustal thickness or (ii) by 
postulating that the strength of the lower mantle increases 
with depth or with time, both of which are possible for 
certain combinations of composition and thermal evolution. 

A small iron-rich core in the Moon cannot be excluded 
by the Moon's mass and moment of inertia. If such a core were 
molten at the time lunar surface rocks acquired remanent 
magnetization, then thermal history models with initially cold 
interiors strongly depleted in radioactive heat sources as a 
primary accretional feature must be excluded. Further, the 

40 

presence of K in an Fe-FeS core could significantly alter 
the thermal evolution and estimated present-day temperatures 
of the deep lunar interior. 



1. INTRODUCTION 

Our knowledge of the chemical and physical properties 
of the Moon is both extensive and detailed. Samples ana- 
lyzed from 11 lunar landing sites to date and information 
from earth-based and lunar-orbital geochimical , optical 
and infrared experiments are being combined to give a 
broad understanding of the present state of a large frac- 
tion of the Moon's surface area. The properties of the bulk 
of the Moon and its interior are less well-defined, but by 
no means obscure. By necessity, study of the lunar interior 
must be largely indirect and the conclusions reached can 
rarely be definitive. Nonetheless, the quantity of data 
that can now be brought to bear on the particulars of the 
deep lunar interior is large. In this paper we shall apply 
these data to the questions of internal lunar structure, 
of compositional models for the subsurface lunar material, 
and of possible evolutionary routes by which the Moon might 
have reached its present state. 

For several reasons, the Moon is a particularly suitable 
and perhaps unique, body for investigating present-day consti 
tution and models of evolution, both of which bear on the 
early history of planet-sized objects in the solar system. 
First, the surface to volume ratio for the Moon is large 



relative to that for other terrestrial planets, so that 
observations made at or near the lunar surface have impli- 
cations for a relatively large segment of the moon. The 
lunar crust, for instance, constitutes approximately 10 
percent of the mass of the Moon; the corresponding figure 
for the Earth is about half a percent. Second, the record 
of the early history of the Moon's surface, particularly 
between 3 and 4 billion years ago, is well preserved. 

Finally, because of the Moon's size, the pressures and 
temperatures of the lunar interior are accessible in the 
laboratory , and thus measurements of important physical 
parameters at the conditions throughout the moon can be 
or have been made. 

We begin with a critical summary of the available con- 
straints on lunar constitution and evolution. The foremost 
constraint on the composition of the Moon comes from the 
chemistry of the lunar surface. Downward extrapolation can 
be made on the basis of the seismic velocity structure, models 
of petrogenesis , and considerations of the density of the 
lunar interior. The thermal evolution of the Moon, and the 
related history of melting and differentiation, must be con- 
sistent with the observed episodes of igneous activity at 
the lunar surface, with estimates of the stress history of 



the lunar lithosphere, and with the present-day temperature 
profile inferred from heat flow and electrical conductivity 
measurements. An interesting but unconfirmed additional 
constraint would be the presence of an ancient lunar magnetic 
field of internal origin. 

We next reconsider in detail possible models for the 
thermal evolution in the Moon. The single heat-flow measure- 
ment, if representative of the Moon as a whole* requires the 
average Moon to be more radioactive than presumed in most 
earlier models. Early upwards segregation of the heat sources 
in the Moon, either during a melting episode which followed 
a rapid accretion or as a primary accretional feature, is 
still required. Differences in present-day temperature or 
in strength of the deep lunar interior are perhaps the only 
means by which thermal history calculations can resolve 
whether the high values of radioactivity at the lunar sur- 
face are due to differentiation in an originally homogeneous 
moon or are partly the result of an inhomogeneous accretion 
process . 

Using several estimates of present day temperature and 
the confirmed thickness of a low-density crust, we examine 
density models for the lunar interior that match the Moon's 
mean density and moment of inertia. Models with chemically 
uniform mantles, though marginally permitted by the present 



uncertainty in the moment of inertia of the Moon, are not 
favored. To fit the moment of inertia, the upper mantle 
must be denser than the lower mantle. The magnitude of the 
required density reversal is probably unstable to solid- 
state creep, though there is less problem explaining the 
stress differences if the upper mantle is richer in Ca and 
Al than the lower or if the lower mantle was relatively 
depleted in radioactive elements at the time of a cold 
accretion. If the moment of inertia value holds up, then 
the lower mantle cannot have the same major element chemistry 
as the bulk Moon. 

Finally, we briefly consider a few of the conse- 
quences of an iron-rich core in the Moon. Even if the lunar 
mantle is uniform, a small core is permitted by the moment 
of inertia. If a convecting core were an early feature in 
the Moon, several further constraints on the initial temper- 
atures and heat sources in the Moon would be imposed. 



2. CONSTRAINTS ON LUNAR STRUCTURE AND EVOLUTION 
2.1 Petrology, Chemistry 

The most direct evidence for the internal composition 
Of the Moon comes from the chemistry, mineralogy, and in- 
ferred petrogenesis of the material returned from the lunar 
surface. Analysis of lunar rock and soil samples , together 
with the results of the orbital x-ray and gamma-ray spec- 
trometer experiments, have shown the lunar surface to be 
comprised principally of three rock types: (1) iron-rich 

mare basalts, (2) nonmare basalts , also labeled KREEP 
and norite, rich in aluminum, radioactive elements and re- 
fractory trace elements, and (3) plagioclase-rich rocks 
(WOOD et al . , 1971; HUBBARD and GAST, 1971; LSAPT, 1972). 

Mare basalt samples have been returned from four dif- 
ferent near-side mare sites to date. Experimental high- 
pressure, high-temperature studies (RINGWOOD and ESSENE, 
1970; GREEN et al . , 1971) indicate that mare basalts pro- 
bably formed by variable amounts of partial melting at 
depths between 100 km and 500 km in a lunar mantle com- 
posed of pyroxene and possibly olivine. The sugesstion 
(GAST et al . , 1970; HASKIN et al . , 1970; PHILPOTTS and 
SCHNETZLER, 1970), based on the abundances of rare-earth 



elements and on the negative europium anomaly in particular, 
that mare basalts were derived from a source rock rich in 
plagioclase is in apparent conflict with the high-pressure 
melting experiments mentioned above. Alternative explana- 
tions, not involving plagioclase in the source rock, for 
the europium anomaly have been advanced by GOLES et al . 

(1970) and RINGWOOD and ESSENE (1970) , and experimentally 
supported by GRAHAM and RINGWOOD (1971). 

Several modes of origin for nonmare basalts have been 
proposed, including a small degree of partial melting of 
a shallow, plagioclase-rich source (HUBBARD and GAST, 1971) , 
more extensive melting of a plagioclase-pyroxine-olivine 
crust followed by crystal fractionation (WALKER et al . , 

1972) , and impact melting of the upper several kilometers 
of a plagioclase-rich crust (GREEN et al . , 1972). GREEN 
et al . (1972) have argued on chemical and mechanical grounds 
that the first two explanations are improbable. A similar 
thread in each of the above hypotheses is the notion of 
abundant plagioclase in the near- surface regions of the 
Moon at the time that nonmare basalts were formed. Nonmare 
or KREEP basalts may constitute a relatively minor, though 
important, fraction of lunar surface material. The gamma- 
ray spectrometer results (ARNOLD et al . , 1972) have so far 
indicated only one major area on the Moon, the Oceanus 
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Procellarum-Mare Imbrium region, where the surface radio- 
activity is sufficiently high to suggest a sizeable KREEP 
component. About 20 percent KREEP and 80 percent mare 
basalt would match the data for that region (LSAPT, 1972) . 

Several authors (WOOD et al . , 1970; ENGELHARDT et al ., 
1970; KING et al . , 1970; SHORT, 1970; SMITH et al . , 1970a) 
hypothesized, on the basis of plagioclase or plagioclase- 
rich fragments in the lunar soil, that plagioclase, pre- 
dominantly anorthite, is an important constituent of the 

• . • f " . ■ ", : , .. 

lunar highlands. This suggestion was confirmed by the Apollo 
15 and 16 x-ray fluorescence experiment (ADLER et al ., 1972a, 
b) which demonstrated significantly greater amounts of alu- 
minum, most likely in feldspar, in the highlands than in 
mare material. Several of the highland regions exhibited 
Al/Si ratios appropriate to anorthositic gabbros and gabbroic 
anorthosites. The anorthosites and plagioclase-rich rocks 
owe their origin to complete or nearly complete melting of 
the outer regions of the Moon (WOOD et al . , 1970, 1971; 
HUBBARD et al . , 1971) . 

2.2 Seismic Velocity 

The lunar seismic velocity distribution, obtained 
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from the travel times and amplitudes of seismic waves from 
artificial impacts on the Moon's surface, can be used both 
to delineate chemical layering and to infer roughly the com- 
position in the outer 100 km of the Moon. The compressional 
velocity profile in the Fra Mauro region of the Oceanus Pro- 
cellarum, constructed from data through the Apollo 16 mission 
(TOKSOZ et al . , 1972a, b), is shown in fig. 1. The very 
rapid rise of velocity with depth in the uppermost 10 km 
can reasonably be ascribed to a self-compaction effect, so 
that the velocity structure suggests three petrologically 
distinct layers: (i) a 25-km thick upper crust; with com- 

pressional velocity (at low porosity) 5 to 6 km/sec; (ii) a 
40-km thick lower crust, with velocity 7 km/sec; and (iii) a 
mantle, with apparent velocity 8 km/sec. A thin (less than 
20 km, say) high-velocity (9 km/sec) layer may also be present 
between the lower crust and normal mantle. The shear wave 

II 

data support the proposed layered structure (TOKSOZ et al . , 
1972b) . 

Comparison of the velocity profile with laboratory 
measurements of compressional velocity in lunar and terres- 
trial rocks at appropriate pressures allows some estimate 
of the most likely rock types for each layer (TOKSOZ et al . , 
1972a, b) . The velocities in the upper crust fall within the 
range of the velocities of mare and nonmare basalts (see 



fig. 1), and thus the upper crust in the Procellarum region 
is prohably similar in composition to basaltic rocks found 
at the surface. The lower crust has a velocity similar to 
those measured in terrestrial gabbros and anorthosites; 
the likelihood of abundant anorthositic gabbros in the lunar 
highlands and the spectral evidence for highland material 
underlying some mare craters (MCCORD et al . / 1972) make 
such a material a strong candidate for the primary consti- 
tuent of this layer. The high velocity in the layer imme- 
diately below the lower crust, if indeed such a layer is a 
widespread feature on the Moon, is greater than that of any 
common basic or ultrabasic rock on earth. A high-pressure 
form of the aluminum rich lower crust has been mentioned as 
a candidate material with the requisite velocity (ANDERSON 
and KOVACH, 1972; TOKSOZ et al . , 1972a, b) . That the average 
velocity of the Moon's mantle, between about 70 and 130 km, 
is about 8 km/sec permits a wide range of rock types, including 
dunites , pyroxenites and exlogites, as candidates for the 
material of the bulk of the lunar mantle. 

2 . 3 Mass, Moment of Inertia 

Perhaps the only firm constraint on the constitution 
of the deep interior of the Moon is that of the limitations 



on the lunar density distribution imposed by the Moon's 
total mass and moment of inertia. We discuss further below 
the importance of these two quantities and the application 
of. physically reasonable density models to the question 
of lunar internal structure. 

Following KAULA (1971), we shall adopt the values 

3 2 . 

GM =4902.7 ±.l km /sec for the product of the gravita- 
tional constant and the Moon's mass, and R = 1737.5 ± .6 km 
for the mean lunar radius. This gives a mean density of 
3.344 ± .004 g/cm"^. 

The precise value of the moment of inertia of the Moon 

about its axis of rotation is more controversial, though it 

is exceedingly important in limiting possible density models 

(SOLOMON and TOKSflz , 1968; RINGWOOD and ESSENE, 1970; 

UREY and MACDONALD, 1971). We shall use KAULA'S (1969) 

2 

value of C/MR = .402 ± .002, which is based on the 

second degree zonal and sectorial harmonic coefficients in 

the lunar gravitational-field solution of MICHAEL et al . , 

(1970) . A more recent and somewhat more complete analysis 

by MICHAEL and BLACKSHEAR (1972) gave a nearly identical 

value. It should be noted that the JPL group (LORELL, 1970; 

LIU and LAING, 1971), using a different analysis technique, 

2 

consistently reports values for C/MR significantly less 
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than .400. We conclude that Kaula's stated uncertainty 
2 

in C/MR may be a bit low. A more accurate and precise 
value must await better far-side and high-latitude tracking 
data. 


2.4 Electrical Conductivity and Temperature 

Several techniques have recently been used for deducing 
the electrical conductivity structure and related properties 
of the lunar interior, including analysis of the Moon's re- 
sponse to transient magnetic field events from simultaneous 
magnetometer measurements on the anti-solar surface and in 
the free-streaming solar wind (DYAL and PARKIN, 1971; DYAL 
et al . , 1972b) , and analysis of the frequency-dependent 
transfer function of sunlit lunar surface to the free-stream 
magnetic field (SONETT et al . , 1971). The inferred electrical 
conductivity profiles were converted to equivalent temper- 
atures (ENGLAND et al . , 1968) to show that temperatures in 
the outer 700 to 800 km of the Moon must be less than 1000°C 
(SONETT et al ., 1971; DYAL and PARKIN, 1971). The more re- 
cent inversion of DYAL et al . (1972b) allows, in addition, 
temperatures in excess of 1000°C by as much as several hun- 
dred degrees at depths greater than 800 km. The temperature 



curve of Sonett et al . has attracted particular attention. 
Conventional thermal history calculations have shown that 
models which realistically account for lunar igneous ac- 
tivity and probable chemical and physical parameters yield 
present-day temperatures significantly higher than the 
Sonett et al . or Dyal-Parkin values (HAYS, 1972; TOKSflz ejt 
al. , 1972c; REYNOLDS et al . , 1972), unless very efficient 
heat transfer by solid-state convection is invoked (TOZER, 
1972) . 

There are several difficulties with the temperature 
estimates from electrical conductivity, however. REIS Z 
et al . (1972) have argued that neglect of the boundary 
condition assymetry imposed by solar wind flow past the 
Moon may affect the conductivity models of Dyal and Parkin 
and particularly of Sonett and others. SILL (1972) has 
pointed out that at the highest frequencies of the data of 
SONETT et al . (1971) , their assumption of uniform source 
field breaks down. A further complication of the sunlit 
surface magnetometer analysis not considered by Sonett 
and others is the compression of a remanent field, largely 
unknown, by the solar wind (SILL, 1972; DYAL et al . , 1972). 
At least as serious an objection to these temperature esti- 
mates is the observation by HOUSLEY and MORIN (1972) and 



DUBA et al . (1972) that the conductivity of the ENGLAND 
et al . (1968) olivine is controlled by ferric ions not 
likely to be important in the reduced environment of the 
Moon. Recalculation of temperature in the Moon from the 
Sonett et al . (1971) conductivity profile, using recent 
data for a ferric-ion-poor terrestrial olivine, gives 
temperatures in excess of 1000 °C below 350 to 600 km (DUBA 
et al . , 1972) . 

2.5 Other indicators of Temperature 

There are several other observations that bear less 

quantitatively on the question of the temperature of the 

lunar interior. The foremost of these is the measured 

2 

value of heat flow, 33 ± 5 erg/cm sec, at the Apollo 15 
landing site (LANGSETH et al . , 1972). Whether this sur- 
prisingly high value is a good indicator of the mean heat 
flow at the lunar surface is an important unanswered ques- 
tion. The high surface radioactivity in the Imbrium- 
Procellarum region (ADLER et al . , 1972) suggests that 
the average lunar heat flow is probably somewhat lower. 

That at least the outer portions of the Moon are rela- 
tively cool at present is implied by the low seismicity 



of the Moon (LATHAM et al . , 1971) and the absence of 
significant tectonic activity. The persistence for some 
3 billion years of the gravity anomalies associated with 
mascons demands that a cool lithosphere several hundred 
kilometers thick has been a permanent feature on the Moon 
since mare formation (UREY, 1968) . The location of moon- 
quake activity at 800 km depth (LATHAM et al . , 1971) and 
the high-frequency seismic waves observed on the lunar 
surface from those events indicate that the Moon can sus- 
tain shear stresses at considerable depth and that little 
or no partial melting is to be expected in the Moon's 
upper 800 km (TOKSfiz et al . , 1972a) • 

2.6 Lunar Chronology 

Age-dating of lunar samples, particularly dating using 
40 39 

Rb-Sr and Ar -Ar techniques (PAPANASTASSIOU et al ., 1970 
PAPANASTASSIOU and WASSERBURG, 1970, 1971a, b, 1972; WASSER- 
BURG and PAPANASTASSIOU, 1971; LUNATIC ASYLUM, 1970; TURNER 
1970, 1970, 1972; TURNER et al ., 1971; HUSAIN et al ., 1972; 
and others), has provided a well-documented outline of 
lunar igneous activity. The time sequence of activity 
places strong constraints on the Moon's early history and, 
to a lesser extent, on its later evolution (PAPANASTASSIOU 



and WASSERBURG, 1971a), 

An approximate outline of events occurring on the lunar 
surface, taken from the works cited above, is given in fig. 
2. The Rb-Sr fractionation requires a major episode of 
melting and differentiation very near the time of lunar 
formation (PAP ANASTAS SIOU and WASSERBURG, 1971a). This is 
a necessary feature of lunar evolution even if some of the 
chemical inhomogeneities in the Moon are the result of 
primary accretional layering (WASSERBURG and PAPANASTASSIOU , 
1971) . Much of the lunar crust may date from this time, 
although intense bombardment, remelting and brecciation may 
have altered most or all of the original rocks of the lunar 
crust. The record of activity between 4.5 and 3.7 billion 
years (b.y.) ago is sketchy and confusing, but some form of 
igneous activity in rocks 15415 at 4.1 b.y., 12013 at 4.0 
b.y., and Apollo 14 rocks at 3.9 b.y. is well established. 
Certainly the major mare basins must have been created at 
some stage in this time interval; the relative ages of the 
mare basin excavations have been reasonably well resolved 
(e.g. MUTCH, 1970), but the absolute ages are controversial. 
The filling of the maria obviously postdated the formation 
of the mare basins. To date, the well-determined ages of 
mare basalts span a relatively narrow time interval between 
3.2 and 3.7 b.y. ago (WASSERBURG and PAPANASTASSIOU, 1971) . 
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The gravity anomalies attributed to mascons probably origi- 
nated during this period. The failure to discover lunar 
igneous rocks younger than about 3.2 b.y. ago is an impor- 
tant constraint on the Moon's thermal history (PAPANASTASSIOU 
and WASSERBURG, 1971a). If melting occurred since that time, 
it must have been very localized in extent or confined to 
depths in excess of several hundred kilometers. 

2.7 Magnetism of Lunar Rocks 

The lunar crust has a relatively large but extremely 
variable magnetization. This has been verified by steady 
magnetic fields measured on the lunar surface (DYAL et al ., 
1970; 1972b) and by magnetic anomalies observed by the 
Apollo 15 orbital magnetometer (COLEMAN et al ., 1972). Fur- 
ther, many returned lunar samples were found to display a 
fairly stable remanent magnetization (STRANGWAY et al . , 1970; 
NAGATA et al . , 1970; HELSLEY , 1970; RUNCORN et al ., 1970; 
and others). On the basis of these observations, which 
appear to require 1000 - Y fields on the lunar surface during 
the emplacement of rocks spanning a time interval of almost 
1 b.y., several workers (RUNCORN et al ., 1970; PEARCE et al . , 
1971, and others) have suggested the fields are internal in 
origin and thus indicate a molten conducting core between 



3.2 and 4.0 b.y. ago. The debate over this conclusion has 
been heated, and many other explanations of the magnetiza- 
tion, including fields associated with the earth, with the 
solar wind, or with shock or thermal effects, have been 
offered. While the matter is not yet settled, the existence 
of a small iron-rich core in the Moon would have implications 
of great importance for models of lunar evolution. Some of 
these implications are explored later in this paper. 


3. LUNAR THERMAL EVOLUTION 

The thermal history of the Moon has been a popular 
subject, particularly since the Apollo missions (HAYS, 1972; 
WOOD, 1972; REYNOLDS et al . , 1972; TOKSfiz et al . , 1972c; 
MCCONNELL and GAST, 1972; TOZER, 1972; HANKS and ANDERSON, 
1972) . The major reason for the proliferation of papers in 
this field is the ever-growing collection of physical and 
chemical data which better define the range of acceptable 
thermal evolution models. We consider in this section all 
of the constraints on lunar evolution outlined above, and 
present thermal history models which closely fit these con- 
straints, The purpose of such calculations can never be to 



arrive at a definitive version of the thermal evolution of 
the moon. Rather, such models are to be used to test and 
to limit proposed features of lunar origin, lunar history 
and present-day lunar constitution. 

This section represents an extension of our earlier 
paper (TOKSOZ et al ., 1972c). The computational technique 
used here is described fully in the earlier work. Briefly, 
a finite-difference scheme is employed to solve the equation 
for conservation of thermal energy in a spherically symmetric 
moon. The effects of melting, differentiation of radioactive 
heat sources, and simulated convection of molten material are 
included. Except where noted below, all physical parameters 
necessary for the calculation are as given in the original 
reference . 

II 

The thermal history models of TOKSOZ et al . (1972c) 
were designed to match the chronology of lunar igneous ac- 
tivity, including the early extensive differentiation of 
the outer portions of the Moon and the extended episode of 
mare filling, and the evidence cited above for a relatively 
cool lunar interior today. The present-day heat flow on 
the lunar surface predicted in the models was about half 
the Apollo 15 value. We consider models below that are 
richer in radioactive elements but are otherwise similar to 



earlier models in order to place limits on the bulk radio- 
activity of the Moon. 

A separate, but related, question is whether the present 
constraints on thermal evolution sufficiently limit success- 
ful models to either favor or reject the hypothesis of pri- 
mary chemical zoning, in particular of radioactive heat 
sources, during the accretion of the Moon (GAST et_al. , 1970; 
and others) . We conclude that the answer is negative at 
present. 

3.1 Evolution Models to Match Heat Flow 

A single measurement of heat flow is a rather tenuous 
thread on which to hang a theory of lunar evolution. Ex- 
perience on earth has taught us the folly of relying on one 
measurement to characterize even a small region. Still, 
the absence of water and significant tectonic activity 
argues somewhat against a great deal of local variability 
of sub-surface temperature. Perhaps more to the point, the 
Apollo 15 measurement constitutes at least half of the direct 
heat flow information we are likely to have for some time. 
That determination must be regarded at present as the best 
available value and possibly in the near future as one of 
two bounds on permissible values. 



In considering thermal evolution models, several 
features of the models we presented earlier (TOKsfiz et al . , 
1972c) are still required. The evidence for large-scale 
fractionation of at least the shallow portions of the MOon 
requires that the outer few hundred kilometers be at near- 
solidus temperatures at the time the Moon was formed (PAPA- 
NASTASSIOU and WASSERBURG, 1971a). This is most readily 
accomplished by rapid accretion (RINGWOOD, 1966; MIZUTANI 
et al . , 1972; TOKsfiz et al . , 1972c); time scales as short 
as 100 to 1000 years must be involved or else the energy 
released by inelastic collisions would have had time to 
radiate back into space. The only forseeable means of 
avoiding the conclusion of rapid accretion is to postulate 
a strong pulse of heat from tidal dissipation, say from 
close approach of the Moon to the earth. While the energy 
released in the Moon by such an event may be of the right 
order of magnitude for melting (SINGER, 1970) , the hypo-: 
thesis is difficult to model quantitatively and to test. 

We continue to start from initial temperature profiles 
derived from a simple model of the accretion process (HANKS 
and ANDERSON, 1969) , and to maintain the scaling of radio- 
active elements K/U = 2000 and Th/U = 4 throughout the 
Moon. For a specified sequence of discrete stages of dif- 
ferentiation of radioactive elements, we may generate a 



one-parameter family of thermal evolution models to assay 
the dependence of the model features on the bulk concentra- 
tion of uranium in the Moon. 

Two such one-parameter families are represented in 
figs. 3 and 4. The initial temperature distributions for 
the models of figs. 3 and 4 are 100-year accretion profiles 
with base temperatures, respectively, of 0°C and 500 °C. 
Differentiation is simulated by upwards concentration at 
discrete times of all radioactive heat sources to the depth 
of deepest melting; after each differentiation the heat 
source concentration is assumed to decay exponentially with 
depth. The skin depth for the decay law is allowed to de- 
crease with time. The skin depth of the last differentiation 
is set so that the present-day radioactivity at the lunar 
surface equals that of the average of soil samples from mare 
sites, about 1.2 ppm (SLIVER, 1970, 1971, 1972), with the 
exception that the skin depth is never taken less than the 
finite-difference grid spacing of 20 km. For the models of 
fig. 3, the last differentiation was 3.7 b.y. ago; for those 
in fig. 4 it was 3.0 b.y. ago. 

Figs. 3 and 4 give the heat flow and two somewhat ar- 
bitrary measures of evolution and present-day state. The 
lithosphere, as defined in the figures, extends to the depth 
at which the temperature approaches to within 100°C of the 
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solidus. The "primordial core" contains undifferentiated 
material; i.e., at the time of last differentiation the 
primordial core was solid. 

2 

The Apollo 15 heat flow value of 33 ± 5 erg/cm -sec 

(LANGSETH et al . , 1972) places rather narrow limits on bulk 
radioactivity of the Moon. For these famililies of evolution 
models, the average uranium concentration in the Moon must 
lie in the range 50 to 80 ppm. These values, and partic- 
ularly the lower limit, are less than figures cited else- 
where (HANKS and ANDERSON, 1972) , principally because of 
the more complete accounting here of heat of fusion and 
energy transfer by convection during melting. In particu- 
lar, our earlier conclusion that the Moon appears to be 
less radioactive than most eucritic meteorites (CLARK et 
al. , 1967; RIEDER and WANKE, 1969) is still valid. 

Somewhat higher values of bulk radioactivity are re- 
quired to reach a given value of heat flow in the models 
of fig. 3 than those of fig. 4. This is due to the more 
complete differentiation in the latter models. Because of 
the higher initial temperatures and the longer allowed dif- 
ferentiation period, the Moon is completely differentiated 
(i.e,, a primordial core is absent) for all models in which 
the present uranium concentration exceeds about 45 ppb, 



whereas the models of fig. 3 possess a primordial core for 
the entire range of radioactivity considered. A feature of 
all models matching the Apollo 15 heat flow is a relatively 
fixed vlaue for the present lithosphere thickness of 400 
to 500 km. 

The thermal history of a particular member of the family 

of fig. 3, the model with U = 70 ppb and a present-day heat 

2 

flow of 33 erg/cm -sec, is given in Fig. 5. The time sequence 
of melting in the model is shown in fig. 6. It may be seen 
that the model satisfies the major constraints upon lunar 
evolution. The zone of melting deepens with time. The 
lithosphere during the first 2 b.y. thickens at the rate of 
about 160 km/b.y. In particular, the shallowest melting 
progresses from 155 to 245 km depth during the period of 
mare filling, in agreement with the depth of origin of mare 
basalts and with the need for a reasonably thick lithosphere 
to sustain the stresses associated with mascon gravity ano- 
malies. Between 1 and 2 b.y. in the model, a convecting core 
develops due to the undifferentiated and relatively radio- 
active primordial material near the Moon's center. The con- 
vecting core grows to a radius of 1200 km at 3 b.y. and then 
slowly shrinks as the concentration of heat sources decays 
with time. The Moon is molten below 580 km depth at present 
in the model. The feature of a molten convecting core is 



absent in models (such as those of fig. 4) in which the 
Moon differentiates completely and concentrates all of the 
radioactive heat sources toward the surface. 

3.2 Models with Primary Accretional Zoning 

The major features of the thermal evolution models 
discussed above can also be satisfied by models in which 
the decreasing concentration of radioactivity in the Moon 
is a primary feature associated with a chemically inhomo- 
geneous accretion. Shown in figs. 7 and 8 is a thermal 
history model following closely one proposed by MCCONNELL 
and GAST (1972). The radioactivity is stratified: if no 

differentiation had been allowed, the present-day abundance 
of uranium would be 120 ppb in the uppermost 200 km, 50 ppb 
in the next 250 km, and 10 ppb in the remainder of the Moon. 
The average present-day uranium concentration is thus 53 ppb 
The ratios K/U and Th/U are held constant. The initial 
temperature profile is from a Hanks-Anderson accretion model 
with a base temperature of 0°C and an accretion time of 
1000 years; the temperature barely exceeds the solidus be- 
tween 40 and 140 km depth. Differentiation of radioactive 
elements is accomplished in four discrete steps between 
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0 and .9 b.y. 

The temperature evolution in the upper few hundred 

kilometers of the moon is similar to models of comparable 

bulk radioactivity discussed earlier. The present-day 

2 

heat flow is 32 erg/cm sec and the Moon is everywhere at 
least 200°C below the solidus. 

The difference in the behavior of the model of fig. 7 
from those discussed earlier is that the center of the Moon 
in this model is tacitly assumed to have been relatively 
depleted in radioactive elements at cold temperatures . 

Thus the lunar center heats slowly and the Moon is every- 
where solid at present. Whether a cold accretion can 
accompany primary chemical layering is an important question 
that must be answered by those advocating such a model. In 
particular such an origin is incompatible with a molten iron- 
rich core at an early stage in lunar history. If questions 
such as these can be resolved and if the Apollo 15 heat flow 
value holds up, then the temperatures in the deep lunar in- 
terior, if discernible, offer the only possible means by 
which thermal history calculations can distinguish primary 
accretional zoning from extensive irv situ differentiation. 

We emphasize that these two alternatives for early 
lunar conditions do not produce very different thermal 
regimes in the present-day Moon. Whatever the mechanism. 



the high, values of heat flow and surface radioactivity 
together with the requirement that the outer portions of 
the Moon are cool and have been so for 3 billion years 
require thermal models in which the radioactive elements 
are strongly enriched toward the lunar surface early in 
lunar history. The extensive Rb-Sr fractionation at a 
similar time requires some near-surface melting (WASSERBURG 
and PAPANASTASSIOU, 1971) , and as yet the extent in depth 
of the early melting episode has not been established. We 
see no compelling reason why a major fraction of the lunar- 
mass could not have participated in this early differentia- 
tion episode. Even a sub-solidus lunar interior does not 
rule out an originally homogeneous Moon. If the Moon is 
initially molten and all heat sources are removed from 
the lunar interior, or if initial temperature and bulk 
radioactivity are high enough so that the entire Moon is 
melted or partially melted within the first 1.5 b.y. and 
if convective heat transport and segregation of radioactive 
elements is permitted in the molten zone after melting is 
only one fourth complete, thermal evolution models satis- 
fying the heat flow can be found which are everywhere solid 
today. (These models have their own difficulties: such 

efficient differentiation of heat sources seems improbable, 
and for some models the near surface temperatures during 



153 


the first 2 b.y. are probably too high to have sustained 
mascon-generated stresses.) Alternatively, convection by 
solid-state creep will give relatively cool present-day 
temperatures in the lunar interior without altering signi- 
ficantly the temperature distribution in the upper few hundred 
kilometers (TOZER, 1972j TOKSOZ et al. , 1972c). 


4. DENSITY MODELS FOR THE LUNAR INTERIOR 

Because the Moon's moment of inertia is very nearly that 
of a homogeneous sphere, little attention has been given 
recently to the density structure of the Moon's deep inter- 
ior. This is understandable, since a first-order approxima- 
tion of uniform density is reasonably accurate for most 
calculations of lunar properties. On the other hand, the 
physical and chemical evidence discussed above, especially 
the mean density and moment of inertia and the particulars 
of the lunar crust, sufficiently limit the range of possible 
lunar density profiles so that several interesting deductions 
may be made about the lunar mantle. 

Below we formulate simple but plausible density models 
that are consistent with the latest values for the Moon's 



mass and moment of inertia, with several estimates of tem- 
perature in the lunar interior, and with the thickness and 
low densities estimated for the lunar crust. It will be 
assumed throughout that density in the Moon is a function 
only of radius. Thus such interesting lateral variations 
in density associated with the departure of the center of 
mass from the center of figure (e.g. , KAULA et al . , 1972), 
with the differences in the Moon' s principal moments of 
inertia (e.g., MICHAEL and BLACKSHEAR, 1972), and with the 
mascons (MULLER and SJOGREN, 1968) will be ignored. 

4-1 Density of Crustal Rocks 

An important aspect of a realistic density model for 
the Moon, particularly when calculating the moment of in- 
ertia, is the inclusion of the relatively thick lunar crust 
indicated by the seismic data. Knowledge of the major rock 
types of the crust and measurements or estimates of the 
density of such materials allows the mean crustal density 
to be inferred with reasonable accuracy. 

The measured densities of mare basalts, corrected for 

3 

porosity, are 3.3 to 3.4 g/cm at standard temperature and 
pressure (KANAMORI et al., 1970; STEPHENS and LILLEY, 1970) . 
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But the maria occupy no more than about one— fifth of the 

Moon’s surface area and the density of the highlands must 

be somewhat less (O'KEEFE, 1968). The density of plagio- 

clase-rich rocks is dominated by anorthite, with density 
3 

p =2.76 g/cm . Measured and model values of anorthosites 

3 

and nonmare basalts are, respectively, 2.8 to 2.9 g/cm 

and about 3.0 g/cm'* (WOOD et al ., 1970, 1971). Assuming a 

four to one ratio of plagioclase-rich highland rocks to mare 

basalts, the mean density of the upper crust of the Moon is 

about 3.0 g/cm 3 (cf. TURXEVICH, 1971). We shall adopt this 

value in calculations below. 

The density of the lower crust is less certain, and we 

shall treat this quantity as a free parameter in the range 

■ 3 

2.8 (nearly pure anorthite) to 3.4 g/cm (mare basalt). We 

favor a narrower range of uncertainty, perhaps 3.0 to 3.1 
: 3 

g/cm , for several reasons. The terrestrial gabbros, norites 

and anorthosites with compress ional velocities matching the 

velocity of the lower crust have densities in the ranges, 

3 3 

respectively, 2.89 to 3.05 g/cm , 2.93 to 3.06 g/cm , and 
2.71 to 2.81 g/cm 3 (PRESS, 1966; ANDERSON and LIEBERMANN , 
1966). Further, the presence of a lower crust suggests that 
mascon gravity anomalies are more likely associated with 
relief along the interface between mare basalts and a less 
dense lower crust (CONEL and HOLSTRUM, 1968) than with relief 
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along a lunar "Moho" (WISE and YATES, 1970) .; Thus a lower 

3 i 

crustal density of 3.0 to 3.1 g/cm is quite reasonable; 
the higher value might be weighted slightly ;by the require- 
ment that the lower crust be more dense than mare basalt 
liquids at near-liquidus temperatures {see WOOD et al . , 1971). 

4 . 2 Temperature '» 

Density is a function of both pressure and temperature. 

The latter variable is especially important in a body as 
small as the Moon, and yet as the discus sioh above has indi- 
cated, temperature is only loosely constrained by several 

■ , • ’• ■ s •' : : 

r ■ 

various observations , many of them qualitative in nature. 

We treat the temperature distribution as imprecisely 
known, and calculate density models for several temperature 
profiles. Four such profiles are shown in fig. 9. Profile 
A is the present-day temperature distribution from fig. 5; 
the surface heat flow matches the Apollo 15 value. Profiles 
B and C are two present-day temperature distributions from 
TOKsfiz et al . (1972c); the total range of final temperature 
profiles in the Toksflz et al . thermal history models is indi- 
cated by the shaded region if fig. 9 (only models with bulk 
uranium concentration equal to 23 ppb are included). Profile 
B (fig. 5 of Toksflz et al .) and profile C (fig. 12 of Toksdz 



et al .) have the same initial temperature curve and bulk 
radioactivity; they differ in that the thermal evolution 
model resulting in profile C is based on the assumed 
dominance of solid-state convection as an energy transport 
mechanism below the lithosphere. Profile D is the temper- 
ature model of SONETT et al . (1971), roughly extrapolated 
to the Moon's center, and is included primarily for his- 
torical reasons as a representative "cool'' Moon. Whatever 
life the model has is due to the vigorous arguments of 
TOZER (1972) that efficient energy transfer by solid-state 
creep would produce a present-day temperature profile not 
very different from that of Sonett and co-workers. 

4.3 Construction of Density Models 

A number of different density models have been calcu- 
lated for layered, compressible, warm Moons. For each model, 
the equations of hydrostatic equilibrium and conservation of 
mass are integrated numerically using a fourth order Runge- 
Kutta method. The integration step-size in each layer is 
chosen so that the total mass and moment of inertia of the 
Moon can be calculated to four-figure accuracy. Temperature 
is taken from an assumed profile. Densities at each discrete 
depth are corrected for compression and thermal expansion. 
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No explicit account is made in the density calculation of 
partial melting where present in a temperature model. 

In the lunar crust, compressibility and thermal expan- 
sion are assumed constant and are taken from reported meas-r 
urements, corrected to low porosity, on mare basalts (STEPHENS 
and LILLEY, 1970, 1971; BALDRIDGE and SIMMONS, 1971) . Values 
used are, respectively, 1.4 Mbar~* and 2.2 x 10 _5 °C~^. In 
the Moon's mantle, the pressure and temperature dependence 
of compressibility and thermal expansion are included. All 
parameters used are those of forsterite (SKINNER, 1962; 

KUMAZAWA and ANDERSON, 1969). 


4 * 4 M odels wit h^moge neous , jlantles_ and min eralogically 

Suppose that the lunar mantle is chemical Inhomogeneous , 
at least to the extent that the entire Moon from a depth of 


65 km to the center has a uniform density at standard tem- 
perature and pressure. (As used below, "standard" tempera- 


ture is -20 °C, the near-surface temperature of the Moon 


according to LANGSETH et al . , 1972.) If a temperature distri- 
bution in the Moon and densities, respectively and p^, 

at standard conditions for the upper and lower crust are 


assumed, the^thp density p^ of the homogeneous mantle 
follows from the total mass of the Moon, and the moment of 



inertia may be readily calculated. However, since we regard 

the density of the lower crust as not well determined, it is 

2 V 

more proper to state that P 3 and C/MR are functions of 

This is illustrated in fig. 10 for the four different 
temperature models of fig. 9. Also shown is KAULA'S (1969) 
value for C/MR , with the uncertainty indicated by the zone 
of shading. 

Several of the density distributions represented in fig. . 
10 are shown in fig. 11. For profile D, the effects on 
density of temperature and pressure nearly balance out in 
the lunar interior; the total range of density is only about 
1 percent in the mantle. For the warmer temperature profiles, 
temperature dominates pressure in the lithosphere while the 
pressure effect is slightly greater than that of temperature 
in the zone of melting and/or convection. The total variation 
of density due to pressure and temperature in the lunar mantle 
for these temperature profiles is 2 to 3 percent. The values 
of bulk modulus and thermal expansion used for profile B 
are also shown. 

For temperature profiles as cool as profile D (SONETT 
et al . , 1971), fig. 10 illustrates that the Moon's moment of 
inertia cannot be matched by a density model that has a uniform 
mantle and yet satisfies the Moon's total mass. For warmer 
temperature models, the moment of inertia can be fit to 
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within the present uncertainty over limited ranges of 

density values for the mantle and lower crust. For profile 

2 

C, a value for C/MR of .400 is obtained if p 2 exceeds 
3.38 g/cm"*, a condition we previously deemed unlikely. This 
would require the density of the mantle, at zero pressure and 

3 

near-surface temperature, to be about 3:40 g/cm . It is of 

note that this value is very close to the density of the 

3 

model spinel pyroxenite (p = 3.42 g/cm ) proposed by RINGWOOD 

and ESSENE (1970) as the material of the lunar mantle. For 

profiles B and A density models With p 2 greater than 3.18 
3 

and 3.27 g/cm , respectively, and uniform mantles give 

moments within the lower range of acceptable values (KAULA, 

1969). These models have mantle densities, at standard con- 

3 

ditions, in the range 3.42 to 3.45 g/cm . Such a density 
could be matched, for instance, by a Ringwood-Essene pyrox- 
enite plus olivine with 100 Mg/ (Mg + Fe) less than about 
83 (cf. GREEN et al . , 1971). The seismic velocity of the 
pyroxenite or an olivine with greater than 17 mole percent 
fayalite is too low to account for the possible 9 km/sec 
layer (ANDERSON and KOVACH, 1972) , but is allowed by the 
most recent results indicating an 8 km/sec mantle. 

2 

That temperature profile A gives lower values for C/MR 
but higher mantle densities than does profile B is due to 
the fact that higher temperatures in the former profile are 
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confined to the uppermost 600 km. The effect of this tem- 
perature difference may be seen in fig. 11. For the same 
values of p^ and p 2 , the mantle above 600 km is less 
dense and the lower mantle more dense using profile A than 

if temperature curve B is adopted. This illustrates what 

2 

should have been obvious : we cannot make C/MR arbitrarily 

large by increasing the hypothesized temperatures within 

the Moon. The temperature at depth will be limited by the 

appropriate melting curves and/or convection stabilization 

temperatures. A temperature distribution such as profile 

2 

B maximizes the value of C/MR for other density constraints 
fixed. We expect that all temperature curves designed to 

2 

fit the Apollo 15 heat flow, for instance, will give C/MR 
less than that of profile B unless a solidus curve substan- 
tially higher in temperature than the one used here is 
proposed. 

For the temperature profile considered above, none of 

the density models with uniform mantles gives a value of 
2 

C/MR equal to .402. Given the present uncertainty in the 

2 

Moon's moment of inertia, models with C/MR as low as .3995 

must be considered acceptable. For profiles C and D and 

other models invoking solid-state convection (TOZER, 1972) , 

2 

allowable values of C/MR are reached only if the lower 
crust has an unusually high density. Similarly, density 



models using warmer present-day temperatures including 
those consistent with the Apollo 15 heat-flow value, give 
acceptable moments of inertia for densities of the lower 
crust somewhat larger than may reasonably be expected. 

Thus, we consider below density models in which the mantle 
is inhomogeneous. Specifically, we treat two-layered mantle 
models that fit both the Moon's mass and moment of inertia. 
There are plausible physical explanations of such inhomo- 
geneous mantles that bear strongly on the Moon's early 
evolution. 

The conclusion that the mantle of the Moon is inhomo- 
geneous, we should note, is dependent on the values assumed 
for the moment of inertia and the crustal thickness. From 
fig. 10, for instance, it may be seen that a value of 
.398 or .399 for C/MR would favor a mantle uniform in 
density and major element chemistry for most temperature 
models proposed for the Moon. Alternatively, if the crust 
in the Fra Mauro region of Oceanus Procellarum is unusually 
thick, say because of the outpouring of a large volume of 

basaltic lava and the subsequent downbuckling of the crust, 

2 

then a homogeneous mantle might be permitted even if C/MR 
exceeds .400. For instance, models with a 30-km thick 

3 

crust of density (at standard conditions) 2.9 to 3.1 g/cm 
and a homogeneous mantle of density such as to fit the 
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total mass of the Moon give C/MR = .400 to 401 if tem- 
perature profile B is used. j 

* 

* . 

’ • 

4.5 Models with Inhomogeneous Mantles ; 

As a simple first approximation to an inhomogeneous 
lunar mantle , we consider the mantle to be composed of 
two homogeneous layers. If the densities of; the crust and 

t • 

* 

a temperature profile are assumed, then the Densities of 

■ ' ' ■» 
the upper and lower portions of the mantle are single- 

valued functions of the depth to the boundary between the 

; 2 

two layers. To satisfy the mean lunar density and C/MR = 

■ • * 

.402, the upper mantle must be denser than the lower mantle 

at similar pressure and temperature conditions. The precise 

density contrast between upper and lower mantle is a strong 

function of other parameters used in constructing the density 

model. ; 

This is illustrated in fig. 12, where we plot p 3 and 

i 

P 4 , respectively the densities at standard iconditions of 
the upper and lower mantle , versus the normalized radius of 

the lower mantle. Curves are shown for the .several temper- 

• ’ .■ ■ ■ 

ature distributions discussed earlier; densities of the 

■ ' t'.. ■ 

upper mantle generally exceed 3.5 g/cm^; lower mantle 



densities are consistently less than 3.5 g/cm^. These 
compare with a mean lunar density, at standard conditions, 
of 3.35 to 3.42 for the temperature models considered. 
Minimum density contrast occurs at a lower mantle radius 
of .7 to .8 R (350 to 500 km depth) . This minimum contrast 
depends on the temperature model and the assumed lower crust 

3 

density but is in the range .16 to .28 g/cm for the curves 
shown in fig. 10. 

A density reversal in the lunar mantle is, of course, 
gravitationally unstable. Thus for such a reversal to have 
persisted for several billion years requires that the Moon 
maintained considerable strength over that time period, at 
least to a depth somewhat in excess of the upper to lower 
mantle transition. This requirement, in view of the possi- 
bility that the deep lunar interior may convect by solid- 
state creep at moderate temperatures (TURCOTTE and OXBURGH, 
1969; TOZER, 1972), strongly favors a thin upper mantle. 

The limitations on mantle structure may be put in some- 
whate more quantitative form. Imagine that a thin dike of 
lower mantle material penetrates the upper mantle to the 
base of the crust. Then because of the density contrast, 
there will be a difference in pressure between the base 
of the dike and the base of the normal upper mantle. The 
maximum shear stress t at the base of the dike, equal to 



half the pressure difference, is plotted for each temperature 

density profile in fig. 12 as a function of lower mantle 

radius. The smallest values of t are 200 to 300 bars, 

depending on the temperature model; t increases with the 

thickness of the upper mantle, though is still within 10 

percent of its minimum value for R/R^ “ • 9 (upper mantle 

thickness 110 km) . Even the smallest values of these shear 

stresses are 2 to 5 larger than the creep strength of the 

lunar crust inferred from crater dimensions (BALDWIN, 1968) 

and the maximum shear stresses estimated to occur beneath 

mascons (CONEL and HOLSTRUM, 1968; O'KEEFE, 1968; WISE and 

YATES, 1970). Thus the density inversions indicated in 

fig. 12 would be unstable to solid-state creep, a conclusion 

that would be strengthened by considering the additional 

density decrease with depth due to the dominance of thermal 

expansion over compression (see fig. 11). We conclude that 

2 

either (i) a value for C/MR of .402 is untenable and 
KAULA'S (1969) lower limit of .400 is more probably an upper 
bound, (ii) the strength of lunar material increases with 
depth, or (iii) the strength of the lunar lithosphere in- 
creases with time and the density inversion implied by 
C/MR = .402 is a feature more recent than the filling of 
the circular mare basins and than the formation of most 
craters. We consider, these possibilities further below. 



The density contrasts between upper and lower mantles 
in the models of fig. 12 were calculated using C/MR^ - .402. 
If Kaula’s lower limit of .400 is assumed, then the density 
contrast and implied shear stresses will be less for a given 

3 

upper mantle thickness. For instance, if =3.0 g/cm , 

3 

P 2 = 3.1 g/cm and temperature profile B is adopted, a 

3 

thm (20 to 40 km thick) upper mantle with density 3.7 g/cm 

2 

at standard conditions will give C/MR = .400; the lower 

3 

mantle density is 3.43 to 3.44 g/cm for such models. Such 
a thin dense layer fits with the suggestion of a high- 
velocity layer beneath the crust discussed above. Further, 
the maximum shear stresses associated with the gravitational 
instability are only 40 to 80 bars. It is not improbable 
that the upper 100 km of the Moon could have supported 
such stresses for the required several billion years. 

4.6 Possible Explanations of a Density Inversion 

The boundary between upper and lower mantle, if valid, 
represents a chemical and/or raineralogical transition. At 
least three possible explanations, not completely indepen- 
dent, of such a transition may be envisioned: (1) The 

upper mantle consists of a high-pressure version of the 
lower crust (cf. ANDERSON and KOVACH, 1972). The crust- 



mantle boundary is a phase boundary associated with the 
breakdown of plagioclase. (2) The upper mantle is the 
residuum left behind after melting and fractionation of the 
lunar crust. The upper mantle-lower mantle boundary repre- 
sents the lower limit of melting in an early period of lunar 
igneous activity. (3) As a result of a differentiation pro- 
cess involving the bulk of the Moon, the upper mantle is 
richer in Fe (SMITH et al . , 1970b) or in Ca and A1 (in 
high-pressure phases) than the lower mantle. 

The possibility that the upper mantle is merely a 
high-pressure form of the lunar crust is explored in fig. 13. 
Shown there are the stability fields of mineral assemblages 
for mare basalt (RINGWOOD and ESSENE, 1970) and a terres- 
trial gabbroic anorthosite (GREEN, 1970) ; also shown are 
the ranges of reported stability boundaries for the plagio- 
clase end members. All of these curves are extrapolations 
from higher temperatures and pressures and must be considered 
approximate at best. Superimposed on the phase boundary dia- 
grams are several possible temperature curves for the shallow 
portions of the Moon; also indicated are the pressures at 
the base of the upper and lower crust. 

It is immediately apparent that if a high density form 
of the aluminum-rich lower crust constitutes the upper mantle 



then the crust-mantle interface is probably not an equili- 
brium boundary. It might, however, be a rate-process boun- 
dary. Though the transition zone for the plagioclase- 
rich material in fig. 13 is wide, it is not necessary that 
the accompanying density increase be a gradual one uniformly 
distributed over the transition interval- The work of ITO 
and KENNEDY (1970) on the basalt-eclogite transition in an 
olivine tholeiite suggests that the largest density jump 
within the transition of the calcium-rich lunar material 
will occur at pressures near to or somewhat less than the 
anorthite breakdown boundary. Thus the crust-mantle 
boundary could conceivably be a phase boundary with a sub- 
stantial density contrast if a rate-limiting temperature 
curve intersects the P-T curve of the largest density in- 
crease at a pressure corresponding to 65 km depth. The 
metastable boundary would have been frozen in when the 
temperature at 65 km cooled below the rate-limiting value, 
an event which might have occurred relatively recently in 
lunar history. 

If the upper mantle and lower crust do not have the 
same composition, but the upper mantle is richer in A1 and 
Ca than the lower mantle, then the density contrast be- 
tween upper and lower mantle may still be due to the insta- 
bility of feldspar in the lunar interior. To illustrate 



this we have drawn in fig. 5, superimposed on the thermal 
history diagram, the approximate boundary between feldspathic 
and spinel pyroxenite for the mare-basalt parent material 
proposed by RINGWOOD and ESSENE (1970). The deduction to 
be drawn from the figure is that the spinel pyroxenite 
assemblage would have become stable in the upper mantle 
only 1 to 2 b.y. after lunar origin. The earliest phase 
change would have occurred at the base of the lithosphere, 
and it is possible that gravitational sinking of chunks 
of the transformed lithosphere into the underlying molten 
zone maintained the lithosphere thickness at a nearly constant 
value for some time. The transformation to the dense 
assemblage for most of the upper mantle might then have 
occurred later than suggested by fig. 5 and thus the shear 
stresses of fig. 12 might be marginally permissible. : 

If both the crust and upper mantle are richer in A1 
and Ca than the lower mantle, then the Moon must have under- 
gone a complicated differentiation process involving most 
if not all of the lunar mass either during or shortly 
after accretion. Suppose, however, that the higher density 
upper mantle is residual material remaining after the melting 
episode during which the lunar crust was formed. Then a 
natural question to ask is whether the lower mantle has 



the same major element composition as the average of crust 
and upper mantle, i.e. whether an evolutionary sequence 
calling for accretion of the Moon from material homogeneous 
in density and differentiation in situ of the crust and 
upper mantle is permitted by the constraints of mass and 
moment of inertia. A requirement of such an evolution is 
that the density of the lower mantle equal the average 
density of crust and upper mantle and, thus, the mean den- 
sity of the Moon, all densities measured at standard con- 
ditions. The mean density of the Moon at standard conditions 
is indicated by an arrow in the lower left of fig. 12 for 
each temperature profile of fig. 3. The lower mantle density 
P 4 is at least . 7 percent less than the mean lunar density 
p for profiles C and D. For profiles A and B, p^ lies 
within . 5 percent of p if R^/R exceeds, respectively, 

.90 and .86, corresponding to depths to the lower mantle 
of 170 and 240 km; the upper mantle density must exceed 
3 . 6 g/cm^ . 

The calculations above have neglected two phenomena 
which affect the comparison of p^ with p: melting and 

solid-solid phase transformations. If portions of the 
lower mantle are molten, then we should raise p^ before 
comparison with the average density of the presumably 
solid crust and upper mantle. (This latter density will 
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still equal the value of p shown in fig. 12, but the mean 
lunar density at standard conditions will not.) Limited data 
indicate that the change in density of a rock upon melting at 
zero pressure is about 10 percent (SKINNER, 1966); because of 
the larger compressibility of the liquid phase (MU RASE and 
SUZUKI, 1968), this density change will decrease with pressure. 
Estimates of the volume change at melting obtained by combining 
thermodynamic data (ROBIE and WALDBAUM, 1968) with the slopes 
of the melting curves for feldspars, pyroxenes and olivines 
(BOYD and ENGLAND, 1963; BOYD et al . , 1964; AKIMOTO et al . , 

1967) lie in the range 3 to 9 percent at 20 kb pressure. Thus 
at depths below 400 km in the Moon a volume change of perhaps 
no more than 6 percent upon melting should be expected. Further, 
it is unlikely that there are extensive regions of complete 
melting; heat transfer by convective transport would quickly 
cool a thoroughly molten zone. We expect that 20 to 25 percent 
melting by volume is a generous upper limit to the extent of 
melting. For larger melt concentrations, upwards migration 
of the liquid phase would be anticipated. To sum up, it is 
unlikely that melting will affect density by more than about 
1.5 percent; because none of our thermal models has a present- 
day temperature profile that shows melting shallower than about 
500 km, the effect of melting on lower mantle density will be 
to increase our estimate of for profiles A and B by about 
.5 to 1.5 percent, depending on the choice of depth to the 
lower mantle . 

The effect of solid-solid phase changes on p 4 is 
critically dependent on the assumed composition. The basalt- 
eclogite transition in mare basalts is accompanied 



by a density increase estimated to be between 7 and 11 
percent (RINGWOOD and ESSENE, 1970; GREEN et al . , 1971). 

The more plagioclase-rich highland and lower crustal material 
should have a larger density jump. The breakdown of pure 
anorthite involves an increase of density of 27 percent 
GREEN (1970) estimated that a terrestrial gabbroic anortho- 
site should give an 11 percent density rise after disappearance 
of plagioclase. A reasonable estimate for the density dif- 
ference between the high pressure and low pressure forms 
of average lunar crustal material is 10 percent. 

The lunar crust constitutes 10 percent of the mass and 
11 percent of the volume of the Moon. If the lower mantle 
is isochemical with a mix of crust and upper mantle, then 
the correction to p 4 which we must apply due to pressure- 
induced phase changes depends on the thickness of the upper 
mantle. We assume the upper mantle is devoid of Al and Ca; 
then the density increase in the lower mantle due to the 
phase- transformed crustal material is the correction to be 
considered subtracted from p^. Lower mantle radii of . 9R, 

.8R, and i7R would imply, for instance, that lunar crustal 
composition amounts to, respectively, 39, 22 and 16 percent 
by volume of the lower mantle. Thus p 4 should be lowered 
by, respectively, 3.9, 2.2 and 1.6 percent for the above 
values of Rj. 



Let us reconsider whether can equal the average 

density of crust and upper mantle if all densities are 
measured at standard conditions. For temperature profiles 
C and D, there is no correction to due to partial 

melting, but should be lowered if solid-solid phase 

transformations are considered. Thus P^ must be less than 
the average density Of crust and upper mantle by about 4 
percent or more for these temperature models. For profiles 
A and B, corrections due both to partial melting and to 
solid-solid transitions must be made. Following the argu- 
ments given above, we corrected p^ for these effects for 
the density models of fig. 12 associated with temperature 
profile B. This corrected (see fig. 12) shows a 

relative maximum at about R 3 / R * • 7 , but for all values 
of r 3 / r / P 4 is at least 2 percent less than the average 
density of crust and upper mantle. A similar conclusion 
can be made for profile A models. 

Thus for all temperature models considered, a value 
2 

for C/MR of .402 can be met only if the lower mantle 
of the Moon is less dense at standard conditions than the 
combined crust and upper mantle. This could be achieved 
by smaller amounts of A1 and Ca or of Fe in the lower 
mantle than in the average crust-upper mantle material. 

In either case, the composition of the lower mantle does 



not match the bulk composition of the Moon and extensive 
differentiation of the Moon, either during or somewhat 
later than the time of accretion, is required. This con- 
clusion has been reached by other workers (HUBBARD and 
GAST, 1971) using somewhat different reasoning. 

That primary accretional zoning of major elements in 
the Moon is a strong possibility lends some support to 
thermal evolution models with primordial zoning of radio- 
active elements. Further, if the contrast in heat sources 
between a radioactive upper mantle and a lower mantle poor 
in radioactivity is sufficiently large, then temperature 
in the lower mantle might have remained cool throughout lunar 
history. In the model of fig. 7, for instance, the Moon 
below 450 km has never been molten and present-day tem- 
perature decreases with depth throughout much of the lower 
mantle. This type of model may offer a solution to the 
stresses implied at the base of the upper mantle by inho- 
mogeneous density models designed to match the Moon's 
moment of inertia. The strength of the lower mantle in 
such a thermal model has been greater than that of the lunar 
crust for most of the Moon's history and has consistently 
been an increasing function of depth below 400 to 700 km. 
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ON A CORE 

The existence of a predominantly-iron core, and the 
requirement that it be molten from 3 to 4 b.y. ago, would 
place interesting constraints on early lunar history. In 
particular, for the lunar center to be at a temperature 
above, say, the Fe-FeS eutectic (BRETT and BELL, 1969) 
at 4 b.y. ago, either the Moon must have formed at such a 
temperature to begin with or have formed at somewhat cooler 
temperatures (~ 500 °C) and have begun with substantial 
radioactivity at depth. In either case. Moon-wide differ- 
entiation is strongly suggested. 

A thermal history model such as that of fig. 5 is 
consistent with this idea. Temperatures .5 b.y. after 
lunar origin exceed the Fe-FeS eutectic, also shown in the 
figure, throughout most of the Moon. Thus core formation 
by melting, sinking and aggregation of iron-rich droplets 
is achieved prior to the formation of the oldest lunar 
rocks possessing remanent magnetization. Thermal history 
models such as that of fig. 7, with initially cold interiors 
depleted in radioactive elements, are inconsistent with an 
early molten core. A feature common to most, if not all, 
thermal evolution models yielding formation of an Fe-FeS 
core prior to 4 b.y. ago is that the core should still be 



molten or partially molten today. Thus depending on the 
sulfur content of the core, one might have to invoke some 
argument other than core solidification (e.g. decay of 
heat sources) to explain the lack of appreciable dipole 
magnetic field associated with the Moon at present. 

The presence today of even a pure-iron core is not 
precluded by the moment of inertia. Limits can be placed 
on the size of such a core, but the limits depend upon 
several assumed properties. Assuming values for the den- 
sities of the core and the upper and lower crust, and 
adopting a present-day temperature profile, one may calcu- 
late the density of the mantle, assumed homogeneous, and 
moment of inertia from the total mass of the Moon. This is 

done in fig. 14 for two rough limits. In the optimistic 

4 

curve, the lower crust is given a density of 3.4 g/cm , 
the temperature distribution is taken from a thermal evo- 
lution model with initially molten Moon with present-day 
bulk uranium concentration 23 ppb, and the core is taken to 
be 12 percent less dense than iron (modeled after the earth' 
outer core, see PRESS, 1970). In the pessimistic curve, the 

3 

lower crust is given a density of 3.0 g/cm , the temperature 

distribution is that of SONETT et al . (1971) , and the core 

2 

is taken to be pure iron. The predicted C/MR for the 
pessimistic model is too low and requires an inhomogeneous 



mantle even if no core is present. The optimistic model 

allows a small core (R /R = .14, or .005 of the lunar mass) 

even with a homogeneous mantle. For such a core, the 

mantle density (at standard conditions) is lowered about 

.3 percent from the core-free value. Using a temperature 

model with higher radioactivity so as to satisfy the Apollo 

15 heat flow would allow an even smaller core (R /R = .09) 

c 

without postulating an inhomogeneous mantle. 

If the Moon contains a small iron-rich core with 
additional sulfur, then the suggestion of LEWIS (1971) and 
HALL and MURTHY (1971) that potassium would have to be en- 
riched in an Fe-FeS melt relative to silicates has important 
implications for thermal evolution. This is illustrated 
in fig. 15. in fig. 15a is shown the thermal evolution 
for an initially molten moon with U = 11 ppb, K/U = 2000, 
and Th/U = 4 . All radioactive elements are concentrated 
toward the lunar surface shortly after the time of lunar 

formation. The model predicts a solid Moon today and a 

2 

surface heat flow of 10 erg/cm -sec. In fig. 15b is shown 
the thermal evolution for a moon with the same initial 
temperature profile and U and Th abundance, but a 
chondritic K/U ratio of 80 000 (WASSERBURG et al . , 1964). 
The near-surface radioactive heat generation is the same 
as in fig. 15a, but 39/40 of the potassium is assumed to 
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remain in a 340 km radius core. The result, while not 

profound, is spectacular. Heat from the additional potassium 

is sufficient to keep most of the Moon melted throughout 

its history. The Moon is molten below 300 km at present 

2 

and the surface heat flow is 17 erg/cm -sec. These calcu- 
lations are not intended to be realistic: the Moon is not 

chnodritic and an iron core, if present, should probably 
not be currently molten. Rather these results are intended 
to demonstrate the importance of properly accounting for all 
heat sources in performing thermal history calculations. 


6. CONCLUSIONS 

A wide assortment of data now constrain models of the 
composition, structure, and evolution of the Moon's interior. 
We have quantitatively examined some of these constraints 
and the suite of lunar models consistent with them. Our 
findings include the following: 

(1) Models of the thermal evolution of the Moon that fit 
the chronology of surface igneous activity, the need 
for a cool, rigid lithosphere since the formation of 
mascons, and the surface concentrations of radioactive 
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elements all involve extensive differentiation and 
upwards concentration of radioactive heat sources 
early in lunar history. This differentiation may be 
either a primary accretional feature or a product of 
nearly Moon-wide melting. 

2 

(2) If the Apollo 15 heat-flow value of 33 ± 5 erg/cm -sec 
(LANGSETH et al. ( 1972) is representative of the Moon, 
then the average concentration of uranium in the Moon 
is currently 65 ± 15 ppb, assuming the Th/U and K/U 
ratios of lunar surface rocks hold throughout the 
interior. Such models have lithospheres 400 to 500 

km thick at present and may or may not be partially 
melted at greater depths. We caution, however, that 
rigid application of surface heat flow as a constraint 
on the Moon's thermal evolution should await the 
Apollo 17 heat-flow experiment. 

(3) Thermal history calculations are only marginally useful 
for distinguishing between primary accretional zoning 

of radioactive heat sources and an early melting episode 
including extensive differentiation. The major difference 
between the two sets of models is the temperature of 
the innermost 1000 km, a region of the Moon for which 



independent information from seismic velocity or 
electrical conductivity studies is currently lacking. 

(4) Density models for the Moon, including the effects of 
pressure and temperature, can be made to satisfy the 
mass of the Moon, KAULA'S (1969) value for C/MR^ of 
.402 ± .002, and the presence of a lOw-density crust 

some 65 km thick (TOKSOZ et al . , 1972a, b) only if the 
mantle of the Moon is chemically inhomogeneous. Spe- 
cifically, the density of the upper mantle must exceed 
the density at similar pressure and temperature of the 
lower mantle by a minimum of 5 percent. The need for 
an inhomogeneous mantle can be avoided if C/MR is 
less than .400 or, alternatively, if the average 
lunar crust is considerably thinner than in the Fra 
Mauro region of Oceanus Procellarum. For models with 
homogeneous mantles, the mantle density is equal to 
that of a pyroxenite or olivine pyroxenite, in agree- 
ment with seismic velocity (TOKSdz et al . , 1972a, b) 
and petrological (RINGWOOD and ESSENE, 1970) consider- 
ations. For the more likely models with inhomogeneous 

mantles , the upper mantle density probably exceeds 
'3 

3.5 g/cm . If the upper mantle was the source region 
for mare basalts, then considerations of mantle 
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mineralogy based on assumed densities (e.g. GREEN 
et al . , 1971) will need revision. 

(5) The density inversion suggested for the mantle of the 

2 

Moon is gravitationally unstable. To match C/MR = .402, 

the maximum shear stress at the base of a column or dike 

of lower mantle penetrating the upper mantle must be 

2 to 5 times the maximum stresses associated with the 

Moon's gravity field. Such stresses could not generally 

have persisted throughout lunar history, suggesting 
2 

(i) C/MR for the Moon is more properly .400 or less; 
or (ii) the high-density of the upper mantle is a 
feature more recent than the last major episode of 
melting and igneous activity, and is perhaps a result 
of a Ca-Al rich upper mantle passing through a phase 
boundary during cooling; or (iii) the deep interior has 
always been relatively cool and capable of supporting 
reasonably large stress differences, a suggestion which 
can only be satisfied by thermal models invoking pri- 
mary depletion of heat sources in the interior of a 
moon accreted at cold temperatures (GAST and MCCONNELL, 
1972) . 

(6) From the density contrast between upper and lower mantle 



2 

required if C/MR = .402, the lower mantle must be 
less dense at standard temperature and pressure than 
the combined crust and upper mantle. This conclusion 
is strengthened by considering the possible effects 
on density of partial melting and solid-solid phase 
transitions. Thus models of lunar evolution calling 
for homogeneous accretion and in situ differentiation 
of the crust and upper mantle are ruled out, and 
primary or early moon-wide differentiation of the major 
elements (especially Fe, Ca, Al) is required. 

(7) A small iron-rich core and a chemically uniform 
mantle are permitted for some temperature models and 
for some values of crustal density. The radius of 
the core can be no more than 14 percent of the lunar 
radius if the mantle is homogeneous. The presence of 
a heavy core in the Moon exacerbates the problems 

of density inversion and implied stresses required by 
the moment of inertia. 

(8) Confirmation of an iron-rich core in the Moon, and 
the need for it to have been molten between 3 and 4 
b.y. ago to have generated an internal lunar magnetic 
field, W ould have very important implications for 
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lunar thermal evolution. The initial temperature 

in the deep interior of the Moon would have to exceed 

1000 °C or, if cooler, would have to be sufficiently 

rich in radioactive heat sources so that such a 

temperature is reached in the first .5 b.y. Further, 

40 

the likely presence of K in an Fe-FeS core would 
alter the distribution of heat sources commonly 
assumed in thermal history calculations and would 
probably yield a molten core at present. 
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Figure 1. 


Figure 2. 


Figure 3 . 


FIGURE CAPTIONS 

Observed compressional velocity profile in 
the Moon (hatched lines) compared with laboratory 
measurements of velocity as a function of 
pressure in lunar (numbered) and terrestrial 
(named) rocks. The lunar profile is from 
TOKSQZ et al . (1972a,b) > where the sources of 
laboratory data are given in full. 

A schematic summary of those aspects of the 
history of the lunar surface relevant to the 
questions of the constitution and evolution of 
the Moon's interior. Sources for this summary 
are given in the text. Items included near the 
top of the figure are well documented; entries 
are progressively more speculative toward the 
bottom. 

The variation of some features of the Moon with 
average present-day uranium concentration for a 
family of lunar thermal history models . All 
models start from an initial temperature profile 
of TOKSOZ et al . (1972c) and given by a model of 
accretion at 0°C over a 100 year time interval 
with an accretion rate from HANKS and ANDERSON 
(1969) , assume Th/U = 4 and K/U = 2000 throughout 



the Moon today, and include the effects of 
simulated convection and differentiation upon 
melting. The last major differentiation was at 
.9 b.y. after origin. Definitions of primordial 
core and lithosphere are given in the text. The 
box outlines the range of uranium concentrations 
that match the Apollo 15 heat flow of 33 ± 5 
erg/cm 2 -sec (LANGSETH et al . , 1972) . Circles 
denote the specific models which were investigated 
and from which the curves were drawn. 

Figure 4. The variation with average present-day uranium 
concentration of heat flow, and dimensions of 
present-day lithosphere and primordial core for 
a suite of thermal history models starting from 
an initial temperature curve based on accretion 
at 500°C and assuming the last major differentiation 
to have been 1.6 b.y. after lunar origin. All 
other parameters and symbols are identical to 
those in fig. 3 . 

Figure 5. A model of temperature in the Moon as a function 
of time since lunar origin. Time, in billions 
of years, is indicated by the number. adjacent to 
each profile. The initial temperature distribution 
is calculated from a simple model of the accretion 



Figure 6. 


Figure 7 . 


process (TOKSOZ et al . , 1972c); accretion at 
0°C over a 100 year time span are assumed. The 
Moon is partially or completely molten at those 
depths where the temperature profile lies along 
the solidus of mare basalt (RINGWOOD and ESSENE, 
1970). The depth range for complete melting is 
delimited by the small arrows above the solidus. 
Also shown are the Fe-FeS eutectic temperature 
(BRETT and BELL, 1969) and the phase boundary 
between feldspathic and spinel pyroxenite irt the 
model lunar mantle composition of RINGWOOD and 
ESSENE (1970) . 

The extent of melting as a function of time in 
the thermal model of fig. 5. Heavy stippling 
denotes complete fusion; lighter stippling 
indicates partial fusion. 

Thermal evolution in a Moon derived from inhomo- 
geneous accretion. The initial temperature 
profile follows from a presumed cold (0°C) 
accretion over a 1000 year time interval. Primary 
zoning of radioactivity is assumed; see text. 

Other parameters and symbols are identical to 
those of fig. 5. 
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Figure 8. The extent of melting in the thermal model of 
fig. 7. Notation as in fig. 6. 

Figure 9. Some estimates of present-day temperature in the 
Moon. Profile A is the final temperature curve 
from fig. 5. Profiles B and C are from TOKSOZ 
et al . (1972c); the range of present-day 
temperature profiles in their models is indicated 
by the stippled region. Profile D is from 
SONETT et al . (1971). The solidus is that of 
mare basalt (RINGWOOD and ESSENE, 1970) . 

2 

Figure 10. The variation of the moment of inertia C/MR and 
the density p 3 of the mantle as a function of 
the density p 2 of the lower crust for models of 
the lunar density profile with homogeneous 
mantles. The letters A through D correspond to 
the temperature profiles of fig. 9. A density 
of 3.0 g/cm^ is assumed (see text) for the upper 
crust. The shaded region indicates the limits 
of uncertainty in the observed moment of inertia 
of the Moon (KAULA, 1969). 

Figure 11. The variation of density with radius in the Moon 
for several models with homogeneous mantles. 

The letters A through D correspond to the temper- 

3 

ature profiles of fig. 9. A density of 3.0 g/cm 
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Figure 13. 


is assumed for the upper and lower crust. Shown 

also are the values of bulk modulus K (in Mbar) 

-5 -1 

and thermal expansion a (in 10 °C ) used in 
conjunction with the density distribution 
appropriate to temperature profile B. 

The variation of densities and of the 
upper and lower mantle as a function of the 
normalized radius of the lower mantle. Density 
models match the mass and moment of inertia of 
the Moon. Densities of 3.0 and 3.1 g/cm^ are 
adopted for the upper and lower crust, respectively. 
The letters A through D correspond to the 
temperature profiles of fig. 9. Also shown are 
the average density of the Moon at standard 
conditions (arrows at lower left) for each 
temperature distribution and the maximum shear 
stress t implied by the density inversion in the 
lunar mantle (see text). The dotted line is p^ 
for profile B corrected for the effects of 
partial melting and solid-solid phase transformations 
(see text) . 

Several estimates of the temperature distribution 
in the shallow portions of the Moon superimposed 
on potentially important phase transition 
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boundaries. Extrapolations of the pressure- 

* 

temperature boundaries for the first appearance 
(with increasing pressure) of garnet and the 
disappearance of plagioclase are shown as solid 
lines for Apollo 11 basalt (RINGWOOD and ESSENE, 
1970) and as alternating short and long dashes 
for a terrestrial gabbroic anorthosite (GREEN, 
1970). The range of equilibrium boundaries , also 
extrapolated, reported in the literature for 
Albite = Jadeite + Quartz (BIRCH and LE COMTE, 
1960; NEWTON and KENNEDY, 1967; BOETTCHER and 
WYLLIE , 1968) and Anorthite = Grossularite + 
Kyanite + Quartz (HAYS, 1967; HARIYA and KENNEDY, 
1968) are shown as stippled bands. Temperature 
profiles include that of SONETT et al . (1971) 
(dashed) , the range of present-day temperatures 
in the models of TOKSOZ et al . (1972c), and the 
range of temperatures spanned by a number of 
thermal history models generated to match the 
Apollo 15 heat flow (LANGSETH et al . , 1972). 

The effect of an iron-rich core of radius R c on 

2 

the moment of inertia C/MR and the mantle 
density in two models of the Moon with homo- 
geneous mantles. In the "optimistic'' model, 



Figure 15. 


the lower crust is assumed to be relatively 
dense (3.4 g/cm ) and temperature curve B is 
adopted. In the "pessimistic" model, a lower 
crustal density of 3.0 g/cm ^ and temperature 
profile p are used. 
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The effect of K in an Fe^FeS core on thermal 
history. (a) Temperature evolution in an 
initially molten Moon with average U concentration 
= 11 ppb, Th/U = 4 and K/U = 2000 at present. 

(b) Temperature evolution in a Moon similar to 
(a) except that K/U = 80,000 in the bulk Moon, 
but K/U = 2000 everywhere except in a core of 
radius R c - . 2R. 
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